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Abstract 


Because of their high strength and stiffness together with their light weight, fiber 
reinforced composite materials offer great potential for applications particularly in the 
aerospace industry. The weight savings translate into increased performance and 
decreased fuel costs. In addition, the use of these materials avoids the dependence on 
foreign sources for the critical elements needed in the new exotic metallic alloys that 
otherwise might be used. Early research has proven the usefulness of these materials, 
but the need to quantitatively characterize important material properties and develop 
applicable nondestructive evaluation techniques remains. Important physical properties 
(i.e. mechanical, thermal, electrical) need to be measured. Moreover, relationships 
between these physical properties and important engineering properties such as 
strength, residual strength after impact and fatigue loading, and fiber-matrix interfacial 
strength need to be examined to provide a basis for quantitative nondestructive evalua- 
tion of these materials. 

Toward this goal, linear and nonlinear elastic properties have been demonstrated 
to be important physical properties in conventional materials. In particular, nonlinear 
properties are important in the nondestructive determination of applied and residual 
stress (strain) as well as measuring the interatomic bonding forces in crystalline solids. 
Also, several investigations have established a possible relationship between nonlinear 
elastic properties and ultimate strength in aluminum and carbon steel. 

This work presents the theoretical treatment of linear and nonlinear elasticity in 
a unidirectionally fiber reinforced composite as well as measurements for a unidirec- 
tional graphite/epoxy composite (T300/5208). Linear elastic properties were measured 
by both ultrasonic and strain gauge measurements. The nonlinear properties were 
determined by measuring changes in ultrasonic "natural" phase velocity with a pulsed 
phase locked loop interferometer as a function of stress and temperature. These meas- 
urements provide the basis for further investigations into the relationship between 
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nonlinear elastic properties and other important properties such as strength and fiber- 
matrix interfacial strength in graphite/epoxy composites. 


Acknowledgments 


The author would like to thank a number of people without whom this work 
would not have been possible. In the Materials Science and Engineering Department at 
Johns Hopkins University, much gratitude is due Dr. R. E. Green, Jr. As my faculty 
advisor he provided much insight into the technical aspects of this work and was a 
tremendous source of encouragement and inspiration. Dr. Ron Eby and Dr. Jim Wagner 
were of technical assistance and provided invaluable comments on the preparation of 
this essay. Also, the students and staff of this department who provided friendship, 
assistance, and support are deserving of mention and gratitude. 

Numerous people in the Materials Characterization and Instrumentation Section 
(MCIS) at NASA Langley Research Center are also due thanks. Dr. Joseph S. Heyman, 
Head of this Section, provided the initial motivation for this research and aided it to 
completion with generous financial support and technical guidance. Dr. Yost, Sid 
Allison, Dr. Cantrell, Peter Kushnick, Jerry Clendenin and others at MCIS provided 
tremendous assistance by sharing their equipment and expertise to make these measure- 
ments possible. 

Thanks is due to my wife Karen and our families who have provided love, sup- 
port, and understanding throughout this work. 

This research was supported by NASA with a Graduate Student Research Fel- 
lowship Training Grant (NGT 21-001-802) as well as the Johns Hopkins Center for Non- 
destructive Evaluation. Thanks are due Dr. Heyman and Dr. Green respectively in this 
regard. 


v 



TABLE OF CONTENTS 


Abstract ™ 

Acknowledgements v 

List of figures 
List of tables 

I. Introduction 1 

II. Linear Elasticity 20 

A. Introduction 20 

B. Theory 20 

C. Experiment 48 

D. Conclusion 72 

III. Nonlinear Elasticity 75 

A. Introduction 75 

B. Theory 75 

C. Experiment 80 

D. Data Analysis 116 

E. Conclusion H7 

IV. Summary and Conclusions 121 

Appendix A: Linear Elastic Moduli for Symmetries of Transverse 124 

Isotropy and Orthotropy 

Appendix B: Equations for Propagation of Linear Elastic Waves in 129 

Transversely Isotropic and Orthotropic Media 

Appendix C: Analysis of Phase Shifts at the Specimen-Transducer Interface 140 

Appendix D: Computer Program for Pulse Echo Overlap Data Acquisition 145 

and Analysis 

Appendix E: Computer Program for Stress-Strain Data Acquisition and 153 

Analysis 

Appendix F: Nonlinear Elastic Moduli for Symmetries of Transverse 160 

Isotropy and Orthotropy 

Appendix G: Computer Program for Stress Acoustic Constant Data 163 

Acquisition and Analysis 

Appendix H: Computer Program to Determine Least Squares Best Fit 169 

for Nonlinear Elastic Moduli 


References 


184 


LIST OF FIGURES 


Figure 1.1 - (a) Load-Elongation curve for linear elastic 
deformation, (b) Load-Elongation curve for nonlinear elastic 
deformation. 

Figure 1.2 - Load-Time and Elongation-Time curves illustrating 
anelastic deformation. 

Figure 1.3 - Load-Elongation curve for plastic deformation. 

Figure 1.4 - Potential Energy (E) - Separation distance (r) curve 
for interatomic bonding. 

Figure 2.1 - Two dimensional representation of shear strain, solid 
line represents undeformed medium while dashed line is deformed 
medium. 

Figure 2.2 - Diagram of nine stress tensor components. 

Figure 2.3 - Diagram of the rotation of axes 180 degrees around the 
x 3 axis (primed axes are the rotated axes). 

Figure 2.4 - Illustration of variation of stress on a volume element 
of material along the direction. 

Figure 2.5 - (a) Illustration of longitudinal (compressional) elastic 
wave demonstrating directions of particle vibration and propagation, 
(b) Illustration of transverse (shear) wave. 

Figure 2.6 - Diagram of a quasitransverse elastic wave undergoing 
energy flux deviation illustrating the directions of particle dis- 
placements, wave normal, and energy flux vector. 

Figure 2.7 - Illustration of axis designation in a unidirectional 
fiber reinforced composite with respect to the fiber direction and 
the lamina stacking direction. 

Figure 2.8 - Illustration of unidirectional composite with resin rich 
regions between laminae causing the material to have orthotropic 
elastic symmetry. 

Figure 2.9 - Block diagram of pulse echo overlap ultrasonic system. 
Figure 2.10 - Illustration of typical pulse echo pattern. 

Figure 2.11 - Illustration of frequency division of a continuous 



wave (c.w.) source by the Matec model 122B by a factor of five. 

Figure 2.12 - Illustration of an oscilloscope trace of two echoes 
during the pulse echo overlap technique where the echoes are not 
properly overlapped. 

Figure 2.13 - Illustration of two properly overlapped echoes. 

Figure 2.14 - Illustration of the oscilloscope trace of the output of 
a mixer when the input c.w. and pulse sources are of the same frequency. 

Figure 2.15 - Diagram of the compression fixture used during the 
static loading for the measurement of the elastic compliance moduli 
and the uniaxial stress acoustic constants. 

Figure 2.16 - Block diagram of the apparatus used for static measurement 
of the linear elastic compliance moduli. 

Figure 3.1 - Block diagram of the pulsed phase locked loop ultrasonic 
interferometer. 

Figure 3.2 - Block diagram of uniaxial stress acoustic constant 
measurement apparatus. 

Figure 3.3 - Plot of normalized frequency change as a function of 
uniaxial stress for a longitudinal wave propagating along X! 
with stress applied along x 3 . 

Figure 3.4 - AF/F versus a for a longitudinal wave propagating 
along X 2 with stress applied along x 3 . 

Figure 3.5 - AF/F versus <j for a longitudinal wave propagating 
along x 3 with stress applied along . 

Figure 3.6 - AF/F versus a for a longitudinal wave propagating 
along x 3 with stress applied along x 2 . 

Figure 3.7 - AF/F versus cr for a longitudinal wave propagating 
along x 2 with stress applied along Xj . 

Figure 3.8 - AF/F versus a for a longitudinal wave propagating 
along Xj with stress applied along x 2 . 

Figure 3.9 - AF/F versus a for a shear wave propagating along 
x 2 polarized along X! with stress applied along x { . 

Figure 3.10 - AF/F versus a for a shear wave propagating along 



X! polarized along x 2 with stress applied along x 2 . 

Figure 3.11 - AF/F versus a for a shear wave propagating along 
Xj polarized along x 2 with stress applied along x 3 . 

Figure 3.12 - AF/F versus er for a shear wave propagating along 
x 2 polarized along X! with stress applied along x 3 . 

Figure 3.13 - AF/F versus <r for a shear wave propagating along 
x 2 polarized along x 3 with stress applied along x 1 . 

Figure 3.14 - AF/F versus a for a shear wave propagating along 
Xi polarized along x 3 with stress applied along x 2 . 

Figure 3.15 - AF/F versus cr for a shear wave propagating along 
x 2 polarized along x 3 with stress applied along x 3 . 

Figure 3.16 - AF/F versus a for a shear wave propagating along 
X! polarized along x 3 with stress applied along x 3 . 

Figure 3.17 - AF/F versus <r for a shear wave propagating along 
x 3 polarized along x Y with stress applied along Xi . 

Figure 3.18 - AF/F versus a for a shear wave propagating along 
x 3 polarized along x 2 with stress applied along x 2 . 

Figure 3.19 - AF/F versus <j for a shear wave propagating along 
x 3 polarized along x 2 with stress applied along Xj . 

Figure 3.20 - AF/F versus a for a shear wave propagating along 
x 3 polarized along Xj with stress applied along x 2 . 

Figure 3.21 - Block diagram of apparatus used in the hydrostatic 
pressure stress acoustic constant measurements. 

Figure 3.22 - AF/F versus pressure for a longitudinal wave propagating 
along Xi . 

Figure 3.23 - AF/F versus pressure for a longitudinal wave propagating 
along x 2 . 


Figure 3.24 - AF/F versus pressure for a shear wave propagating along 
Xi polarized along x 2 . 

Figure 3.25 - AF/F versus pressure for a shear wave propagating along 
x 2 polarized along Xj . 



Figure 3.26 - AF/F versus pressure for a shear wave propagating along 108 

Xj polarized along x 3 . 

Figure 3.27 - AF/F versus pressure for a shear wave propagating along 108 

x 2 polarized along x 3 . 

Figure 3.28 - AF/F versus pressure for a shear wave propagating along 109 

x 3 polarized along x t . 

Figure 3.29 - AF/F versus pressure for a shear wave propagating along 109 

x 3 polarized along x 2 . 

Figure 3.30 - AF/F versus temperature change for a longitudinal wave 111 

propagating along Xj . 

Figure 3.31 - AF/F versus temperature change for a longitudinal wave 111 

propagating along x 2 . 

Figure 3.32 - AF/F versus temperature change for a shear wave propagating 112 

along x 3 polarized along x x . 

Figure 3.33 - AF/F versus temperature change for a shear wave propagating 112 

along x 3 polarized along x 2 . 

Figure 3.34 - AF/F versus temperature change for a shear wave propagating 113 

along Xj polarized along x 2 . 

Figure 3.35 - AF/F versus temperature change for a shear wave propagating 113 

along x 2 polarized along . 

Figure 3.36 - AF/F versus temperature change for a shear wave propagating 114 

along Xi polarized along x 3 . 

Figure 3.37 - AF/F versus temperature change for a shear wave propagating 114 

along x 2 polarized along x 3 . 

Figure 3.38 - AF/F versus temperature change for a longitudinal wave 115 

propagating along x 3 . 


x 



LIST OF TABLES 

Table 2.1 - Equations for pure mode linear elastic wave propagation 49 

in transversely isotropic and orthotropic media. 

Table 2.2 - Equations for non pure mode, off axis linear elastic wave 50 

propagation in transversely isotropic and orthotropic media. 

Table 2.3 - Sample dimensions and densities for samples used for pure 61 

mode elastic wave propagation. 

Table 2.4 - Sample dimensions for samples used for off axis, non pure 61 

mode wave propagation. 

Table 2.5 - Data for pure mode linear elastic wave propagation. 61 

Table 2.6 - Data for off axis, non pure mode linear elastic wave 64 

propagation. 

Table 2.7 - Measured linear elastic moduli for transversely isotropic 64 

and orthotropic symmetries. 

Table 2.8 - Comparison of measured transversely isotropic stiffness 65 

data with that presented by Kriz and Stinchcomb [45]. 

Table 2.9 - Calculated linear elastic compliance moduli using ultra- 67 

sonic stiffness data for transversely isotropic and orthotropic symmetries. 

Table 2.10 - Measured linear elastic compliance values for transversely 73 

isotropic and orthotropic symmetries. 

Table 2.11 - Comparison of measured compliance values with those 73 

presented by Garber [48]. 

Table 3.1 - Stress Acoustic Constant equations for transversely 81 

isotropic and orthotropic media. 

Table 3.2 - Measured values for the Uniaxial Stress Acoustic Constants. 90 

Table 3.3 - Measured and corrected values for the Hydrostatic Stress 105 

Acoustic Constants. 

Table 3.4 - Measured and corrected values for the Thermal Acoustic 110 

Constants. 

Table 3.5 - Calculated values of third order elastic stiffness moduli 118 

for transversely isotropic and orthotropic symmetries. 


xi 



I. INTRODUCTION 


The theory of nonlinear elasticity and the measurement of nonlinear elastic pro- 
perties of materials have been subjects of great interest for many years. Nonlinear elas- 
ticity is important in many areas of physical measurements including the characteriza- 
tion of interatomic bonding forces and anharmonic behavior of a crystal lattice, nondes- 
tructive measurement of applied and residual stress (strain) using ultrasonic techniques, 
and stress-strain relations in finite strain of materials. This will be discussed in more 
detail later. 

One of the first successful theoretical treatments of the subject of nonlinear elas- 
ticity was presented by Murnaghan [l] in 1937. He extended the theory of infinitesimal 
elasticity to the general case including strains of any magnitude. Birch [2], in 1938, 
pointed out the importance of this theory for the case of seismic waves superimposed on 
the large hydrostatic pressures within the earth. Since that time a number of research- 
ers have presented volumes of material on the subject of nonlinear elasticity. Although 
a detailed review of this information is beyond the scope of this essay, important work 
in several areas will be discussed. First, literature discussing the basic property of non- 
linear elasticity and its effects on ultrasonic wave propagation will be reviewed. This 
will be followed by a discussion of some of the literature involving the development of 
theoretical and measurement techniques as well as some of the applications for non- 
linear elastic research. Also, previous measurements of the nonlinear properties of poly- 
mers and composites will be discussed. 

The study of the mechanical properties of materials is concerned with the defor- 
mation of a material upon application of an external force as well as its behavior upon 
removal of that force. There are three types of mechanical behavior exhibited by 
materials under loading. These are elastic, anelastic and plastic. In elastic deforma- 
tion, the material deforms instantaneously upon application of load and returns to its 
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original shape instantaneously upon removal of load. This behavior can be either linear 
or nonlinear. This is demonstrated in Figs. 1.1 (a) and 1.1 (b) for the simple case of 
uniaxial loading where the elongation serves as a measure of the deformation of the 
material. In anelastic deformation, the material deforms over time upon the application 
of load and returns to its original shape over time upon the removal of load. This is 
demonstrated in Fig. 1.2. When a material undergoes plastic deformation, it deforms 
instantaneously upon application of load but does not return to its original shape upon 
the removal of the load. Fig. 1.3 shows the load-elongation characteristics of plastic 
deformation. 

However, the behavior of real materials can be a combination of these three 
types of mechanical behavior and is generally dependent on the amount of load applied. 
The material may be elastic for small loads but the deformation may be anelastic or 
plastic at higher loads. Also changes in other parameters such as temperature may 
affect the response of the material upon loading. The present work is interested in the 
case of nonlinear elastic behavior where the deformation is recoverable but not linearly 
related to the applied load. 

An understanding of linear and nonlinear elasticity can be gained by considering 
theories of interatomic bonding in molecules such as the Born [3] and Madelung [4] 
theories for ionic crystals. In these theories, the total energy of bonding is the energy 
due to the sum of the attractive and repulsive forces between the atoms multiplied by 
the distance over which they act. The general shape of such a curve is shown in Fig. 1.4 
with the potential energy (E) versus the separation distance (r) plotted. The equili- 
brium separation distance is r 0 where the potential energy is at a minimum Eq . To 
move the atoms closer together or further apart raises the energy of the system and 
thus requires the application of an external force. This force can be determined by tak- 
ing the derivative of E with respect to separation distance for the curve in Fig. 1.4. The 
resulting force versus separation distance curve is then a simple analogy on the atomic 
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scale to the load-displacement curves discussed earlier. 

For small displacements about r 0 , the potential energy curve appears to be and 
can be approximated as a quadratic function. Therefore its derivative is a linear curve. 
That is, the relationship between force and separation distance (load and elongation) is 
linear for small displacements about equilibrium. However, as the amplitude of the dis- 
placement becomes larger, the potential energy curve is no longer quadratic. This 
implies that the force-separation curve is no longer linear which gives rise to nonlinear 
elasticity. 

Another interesting point about Fig. 1.4 is the asymmetry of the potential energy 
curve about r 0 . This provides a simple explanation for the understanding of thermal 
expansion in materials. Since the curve is asymmetric, as the material is heated and 
the atomic lattice vibrations increase, the mean separation distance between atoms 
changes. For the curve in Fig. 1.4, the mean separation distance increases with increas- 
ing amplitude vibrations about r 0 . Thus, in agreement with empirical observations in 
most materials, expansion occurs as the material is heated. 

Because of the small strains (on the order of 1CT* ) imposed on the sample during 
ultrasonic wave propagation, nonlinear elastic behavior is generally neglected in the 
theories of the propagation of ultrasonic waves. However, as pointed out by Green [5], 
there are three cases where nonlinear elasticity may become important in elastic wave 
propagation. First, the amplitude of the wave may become large enough to cause finite 
strains in the material. In the second case, nonlinear behavior may occur when a small 
amplitude wave is superimposed on a large external static stress. Also, nonlinear effects 
may be caused by defects in the material which cause localized regions of finite strain. 

The effects of nonlinear elasticity on ultrasonic wave propagation were also 
pointed out by Green [5]. These include the distortion of finite amplitude sinusoidal 
longitudinal waves. As they propagate, energy from the fundamental frequency is 
transferred into the harmonics that are generated. The effect has been used by 
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numerous investigators to measure the nonlinear properties of solids. Measurements on 
aluminum single crystals as well as NaCI, KC1, LiF and a magnesium aluminum alloy 
were presented by Gedroits and Krasil’Nikov [6] in 1963. Breazeale and Thompson [7] 
reported measurements of harmonic distortion in polycrystalline aluminum also in 1963 
which was followed by numerous other investigations. 

Another effect of nonlinear elasticity on ultrasonic wave propagation is that a 
pure mode nonlinear longitudinal wave may propagate alone, but a nonlinear transverse 
wave must have a longitudinal component. Also, nonlinear transverse waves do not dis- 
tort when propagating in a defect free solid. 


Nonlinear elastic waves can interact with each other in a solid to produce other 
waves. This interaction can also take place with thermal phonons to cause energy loss 
from the elastic wave. The amount of interaction depends on the amplitude of the elas- 
tic waves. Interactions of two nonlinear ultrasonic waves to produce a third ultrasonic 
wave were demonstrated by Rollins [8] in 1963 in fused silica, polycrystalline aluminum 
and polycrystalline magnesium. Rollins, Taylor and Todd [9] demonstrated this 
phenomenon again in 1964 and reported the measurements to be correct to within an 
order of magnitude with predictions based on the measured third order elastic constants 
of polycrystalline magnesium. Since this time, further examinations of nonlinear wave 
interactions have followed. 

Another effect of nonlinearity in elastic wave propagation is that the velocity of 
a small amplitude ultrasonic wave superimposed on a large static stress is dependent on 
the amount of static stress applied. It is this effect which provides the basis for the 
measurement technique used in this research and will thus be discussed in more detail. 
The effect can be accounted for by two simple theoretical explanations. In the first, the 
equation for the velocity of an elastic wave propagating along a long thin rod is recalled 
to be 


v 


2 


(i.i) 
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where v is the ultrasonic phase velocity, 

E is Young's modulus 
and p is the density. 

This is based on the assumption that the material is homogeneous, isotropic and 
behaves in a linear elastic fashion. These assumptions imply that the constitutive equa- 
tion of the material has the form 


<T = Ee 


( 1 . 2 ) 


where is the stress 
and € is the strain. 

If a nonlinear constitutive equation is used such as 

<r = E 1 e + E 2 e 2 , (1.3) 

it follows that the velocity will now depend on strain and therefore stress. 

An alternative approach to understanding this is based on the potential energy 
theory for interatomic bonding discussed earlier. Since the velocity of ultrasonic waves, 
together with the density of the medium give a measure of the modulus (second deriva- 
tive of the potential energy), this modulus can be evaluated as a function of strain (or 
stress) by measuring the velocity as a function of strain (or stress). Also, since for large 
displacements (i.e. finite strain), the potential energy function is not a simple quadratic, 
the modulus and thus the velocity would be expected to change with applied stress. As 
will be seen later, the relationship between velocity and stress should be linear if terms 
up to third order in strain are included in the elastic theory. 

The theoretical development of the equations relating ultrasonic velocity to 
stress began with the original work on finite deformation of materials by Murnaghan [l]. 
This work provided the basis for the theory presented by Hughes and Kelly [10] in 1953 
in which they derived wave speeds as a function of stress in homogeneous, isotropic 


- 10 - 


materials. They considered both uniaxial and hydrostatic compression for different com- 
binations of longitudinal and transverse (shear) mode waves. The following seven equa- 
tions were presented: 


Po v xO = X + 2 \x — —“(61 + 4m + 7X + 10^) (1.3) 

P 0 v y0 = » ~ “ “ n + 3X + M ( L4 ) 

Po v xx = X + 2/i - -£-(21 + x + — — ^“(4m + 4X + 10/z)) ( 1 . 5 ) 

3K 0 /i 

Wh ~ " - i^ (m + ^ + 4X + 4f,) (1 ' 6) 

Po v xy = X + 2/i - -“J-(2l “ + X + 2/i)) (1.7) 

wh = <* - ^ + x + ^ “■ 8) 

^ " - ^ (m - (A i^ )n _ 2X) ,19) 


where \,fi are the second order linear elastic moduli of Lame’ for a homogeneous, 
isotropic material, 

l,m,n are the third order elastic moduli for a homogeneous, isotropic 
material as defined by Murnaghan, 
p Q is the density of the unstressed material, 

K 0 is the bulk modulus, 

P is the uniaxial or hydrostatic compression, 
and v is the velocity of an ultrasonic wave propagating along the x axis. The 
first subscript on v refers to the direction of polarization, while the 
second gives the direction of loading. 0 for the second subscript implies 


hydrostatic loading. 
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They made measurements of the nonlinear constants of polystyrene, Armco Iron and 
Pyrex. Toupin and Bernstein [11] followed this work with an extension of the theory to 
include a material of arbitrary crystal symmetry. This work was published in 1961. 

The next major theoretical developments were presented in a series of papers by 
Brugger [12,13] and Thurston and Brugger [14] in 1964 and 1965. These papers provided 
the formulation and notation for most of the research presented in following years on 
the topic. In the first paper by Brugger [12], a formal thermodynamic definition of the 
higher order elastic constants of a solid was given. The relationships between his 
definition and those used by Murnaghan [l] and other authors were also presented. In 
the paper of Thurston and Brugger [14], the concept of an ultrasonic ’’natural” velocity 
was introduced. The ’’natural” velocity was defined as the initial (unstressed) path 
length divided by the time of flight of the elastic wave. The concept of ’’natural” velo- 
city is of great importance in the measurement of nonlinear properties because it 
requires that the change in only one parameter, the time of flight, be monitored as a 
function of stress in the experiments. In determining normal ultrasonic phase velocity 
changes, the time of flight as well as the propagation distance must be measured mak- 
ing the experimental measurements more difficult and uncertain. Also presented in this 
paper were general equations for the relationship between the ’’natural” velocity 
changes and applied stress for different crystal symmetries, wave modes, and stress con- 
ditions. In the final paper by Brugger [13], explicit equations were presented for these 
relations for a number of crystal symmetries, wave modes and stress conditions. These 
equations allow the calculation of all of the third order elastic coefficients for all of the 
crystal point groups. 

Thus, based on the theory of Thurston and Brugger [14], measurements need only 
be made of the change in "natural” velocity (i.e. changes in time of flight) of ultrasonic 
waves as a function of stress to determine the third order elastic coefficients. However, 
the changes involved are generally quite small and require very accurate ultrasonic 
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techniques. Several such techniques have been developed over the years which make 
these measurements possible. Early measurements by Hughes et al. [15] in 1950 used a 
simple pulse transmission technique. A voltage spike was used to excite a piezoelectric 
ultrasonic transducer which launched an ultrasonic wave into the sample. A receiving 
transducer on the opposite side of the sample converted the elastic wave into an electri- 
cal signal which was displayed on an oscilloscope. The time of flight was measured 
directly on the scope with an uncertainty reported to be +/— 0.03 /is. This led to frac- 
tional uncertainties of larger than one part in 10 3 for samples with a 20 /is travel time. 

Improvements in the ability to measure more accurately the ultrasonic wave 
time of flight followed rapidly. In 1953, McSkimin [16,17] reported on a phase com- 
parison ultrasonic technique. This technique was a pulse-echo technique which used the 
same transducer as a generator and receiver of the elastic waves. A tone-burst (burst of 
several cycles of high frequency sinusoidal waves) was used to excite the transducer and 
thus generate the ultrasonic pulse. The elastic pulse traveled through the sample and 
reflected off of the other side. This occurred a number of times producing an oscillo- 
scope pattern of a number of decaying amplitude pulses which were all similar except 
for the exponential decrease in their amplitude. In McSkimin’s technique he amplified 
and rectified the received echos and displayed them on the oscilloscope. The time of 
flight was determined by varying the frequency of the tone-burst to find frequencies 
where the received echos were all in phase. For this condition the time of flight was 
given by 



where t was the time of flight, 
n was an integer, 

f was the frequency of the tone-burst 
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and 7 was a phase angle, in degrees, correcting for phase shifts at the 
specimen-transducer interface. 

The integer n was determined by measuring two sequential values of the frequency for 
the "in phase" condition. In this case, n was given by 


n 


Af 


(i.n) 


where f x and f 2 were the two measured values of frequency 
and Af was the difference between the two measured values. 

An expression for the phase angle 7 based on standard equivalent electrical transmission 
line theory was presented. It was calculated from known material properties (mechani- 
cal impedances of specimen and transducer and resonant frequency of transducer). The 
overall uncertainty in measuring velocity was reported to be 0.14 %. 

In 1957 Williams and Lamb [18] reported a technique using a pair of ultrasonic 
bursts and a through transmission technique. The first burst was launched through the 
sample and reflected multiple times through the specimen. The second burst was 
launched through the same transducer a short time later so that the first arrival of the 
second pulse at the receiving transducer would coincide with the arrival of the first echo 
of the first pulse. The excitation frequency was then adjusted so that there was cancel- 
lation between the two signals. Analysis was then presented which allowed calculation 
of the travel time in the specimen from the frequency which included corrections for 
phase shifts at the specimen-transducer interface. Analysis of the phase angle also 
included the effect of the coupling medium between the transducer and specimen. The 
measured velocity was considered to be accurate to one part in 10 4 . 

Other techniques developed later such as the pulse superposition method 
(McSkimin [19]) in 1961 and the pulse echo overlap technique (May [20] and Papadakis 
[21]) increased absolute accuracy of velocity measurements under certain conditions to 
several parts per million according to Papadakis [22]. He also reported that changes in 
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time of flight could be measured accurate to parts in 10 7 using these techniques. Details 
of the pulse echo overlap technique will be discussed later. 

The final technique discussed was reported by Heyman [23,24] in 1980. It 
allowed changes in "natural" velocity accurate to one part in 10 7 to be measured. The 
technique uses an ultrasonic instrument known as a pulsed phase locked loop (P2L2) 
which measures changes in resonant frequency of the specimen, transducer, bond compo- 
site resonator system as a function of some other parameter such as stress or tempera- 
ture. This can be shown to be equivalent to measuring changes in "natural" velocity 
(W) as follows. Since 
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where L 0 is the initial length of propagation which is a constant 
and t is the time of flight, 

then 

7 


AW = 


(t ALq 


LoAt) 


(1.13) 


which reduces to 


AW = - 


Lq At 


(1.14) 


since 


ALq = 0. 


Therefore, 


AW _ At 
W t 


(1.15) 


(1.16) 


The normal phase velocity (v) is given by 
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V = 


L 

t 


(1.17) 


where L is the length of propagation which is not a constant. 
Therefore, 

AW = At _ Av AL 

W t v L ' 


However the resonant frequency of the test system is given by 


F 


n 


nv 

2L 


(1.18) 


(1.19) 


where n is an integer denoting the number of the harmonic. 
Therefore, 

AF _ Av AL 
F “ v “ L ’ 

and thus, 

AF = AW = At 
F W t ' 


( 1 . 20 ) 


( 1 . 21 ) 


The way in which the P2L2 measures changes in resonant frequency will be discussed in 
a later section as it is the instrument used for the nonlinear measurements in this work. 

Although most of the work discussed thus far deals with the development of 
theory and measurement techniques for nonlinear behavior of materials, much work in 
the literature points out the importance of nonlinear measurements and their applica- 
tions. One of the first applications of nonlinear elasticity was in the area of materials 
characterization. Nonlinear elastic constants, through their representation in terms of 
interatomic potentials in crystalline solids, provide much information about the nature 
of the bonding of atoms in a crystal. Einspruch and Manning [25] pointed out the 
importance of nonlinear and anharmonic phenomena in solid materials in relation to 
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other important material properties such as thermal coefficients of expansion, Griineisen 
constants and ultrasonic attenuation. They also pointed out the need for a nonlinear 
description of elasticity in theories dealing with regions of high stress or finite strain in 
solids such as around dislocations. 

One of the most useful and important applications of nonlinear elasticity, how- 
ever, has been in the field of nondestructive evaluation of applied and residual stress 
(strain) in materials. This application is based on the fact that ultrasonic wave speeds 
are a function of the stress state of a material due to nonlinear effects. Thus, if the 
relationship between stress and velocity is previously measured for a given material, 
theoretically, the state of stress in the material can be monitored by measuring the 
ultrasonic velocity. A variation of this is called the acoustoelastic effect or acoustical 
birefringence. The acoustoelastic effect arises from the fact that two shear waves, one 
polarized parallel to stress and the other perpendicular to stress, have a difference in 
velocity proportional to the applied stress. Again this is due to nonlinear effects. 
Acoustical birefringence is similar to optical birefringence or the photoelastic effect 
which has been used to determine strain in transparent materials for many years. As 
early as 1959 the acoustoelastic effect was being touted for the measurement of residual 
stress. Benson and Raelson [26] stated that acoustoelasticity could be as effective as 
photoelastic methods with the additional benefit of applications to opaque materials. 

However, in the years since, problems have plagued the application of acoustical 
birefringence for the measurement of residual stress. The effect of velocity difference 
between the shear waves due to material anisotropy is often much larger than changes 
observed due to stress making residual stress measurements difficult. Another problem 
involves energy flux deviation. In an anisotropic body, the energy flux vector (the direc- 
tion of the flow of energy per unit time per unit area) of an elastic wave does not in 
general coincide with the wave normal. This refraction of the wave is different for the 
two shear waves which means the two waves do not follow the same path through the 
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material. It also depends on the degree of anisotropy and may also change as a func- 
tion of stress. These factors make velocity measurements for acoustoelastic determina- 
tion of stress extremely difficult. Also, as pointed out by Hsu [27], although there are 
analogies between the photoelastic effect and the acoustoelastic effect, a comparison of 
the important parameters (i.e. wavelengths, frequencies, velocity changes, etc.) for both 
acoustical and optical wave propagation shows that the acoustic technique is not as 
efficient as the optical technique. 

Nevertheless, numerous investigators are still attempting to overcome these prob- 
lems. Theory has been developed to attempt to account for the effect of special textures 
which has met with limited success in applications. Also, as noted earlier, measurement 
technology in ultrasonics has progressed substantially. The progress in acoustoelasticity 
was reviewed in detail by Pao et al. [28]. The general problem of the experimental 
determination of residual stress using ultrasonic methods has yet to be solved. However, 
the detection of applied stress has been successfully reported by numerous authors 
[24,29,30]. 

Nonlinear elastic material properties may also be important in the nondestruc- 
tive evaluation of important engineering properties of materials. For engineering appli- 
cations, it is important to determine material properties such as strength and residual 
strength after impact or fatigue loading' and thermal cycling. However, most nondes- 
tructive measurements yield only physical properties such as modulus, density or 
coefficient of thermal expansion. To determine the important engineering properties 
relationships must be developed between the measured physical quantities and the 
desired engineering properties. Nonlinear properties may be useful in this respect as 
they are parameters which are significant in large deformations (i.e. near failure strains) 
of materials. To examine this, Heyman et al. [31] studied the relationship between the 
Stress Acoustic Constant (SAC) which is a measure of a mixture of second and third 
order elastic moduli and carbon content in carbon steels. They found the SAC not only 



- 18 - 


seemed to be related to the carbon content but to the strength of the material as well. 
In aluminum, Heyman and Chern [32] examined the relationship between the SAC and 
heat treatment. They found that while the second order moduli were insensitive to the 
different heat treatments used, the SAC was able to differentiate them. Thus, a simple 
test method for detecting improperly heat treated and thus inferior strength aluminum 
was demonstrated based on measurements of nonlinear elastic properties. 

Although previous discussion has considered primarily only metallic and non- 
metallic crystalline solids, numerous studies have presented measurements of nonlinear 
phenomena in polymeric and composite materials. In 1950, Hughes et al. [15] reported 
changes in wave speed as a function of pressure and temperature for polyethylene, 
polystyrene and Lucite. In 1953, Hughes and Kelly [10] presented the three third order 
elastic constants for polystyrene. In 1959, Singh and Nolle [33] measured the velocity as 
a function of hydrostatic pressure and temperature for poly isobutylene. Assay et al. 
[34,35] presented the change in velocity for polymethylmethacrylate as a function of 
stress and temperature. The relationship between stress and velocity became nonlinear 
at high pressures indicating higher order nonlinear effects. Lamberson [36], in 1969, 
reported on the temperature and hydrostatic pressure dependence of velocity in polys- 
tyrene and two composite materials. One of the composites was carbon phenolic while 
the other was a tape wound silica phenolic. He also reported nonlinear stress-velocity 
curves for these materials at high pressures. Zarembo and Shklovskaya [37] reported 
harmonic generation data for polystyrene, plexiglass and rubber in 1971. They also 
examined the effect of hydrostatic and uniaxial stress on harmonic generation in these 
materials. Other authors [38,39] have reported nonlinear measurements on different 
polymeric materials since this time. 

The present work deals with the measurement of the linear and nonlinear elastic 
properties of a graphite/epoxy composite material. More specifically, measurements 
were made on unidirectional laminates of T300/5208 (Thornel 300 graphite fibers in a 
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Narmco 5208 epoxy resin). The choice of a unidirectional lay up was made because it 
offers the highest elastic symmetry and therefore is the least complex orientation to 
analyze. The elastic symmetry of unidirectional fiber reinforced composites is usually 
taken to be that of transverse isotropy if the fibers are randomly distributed in the 
plane perpendicular to the fiber axis. A transversely isotropic material is one which has 
a preferred axis, perpendicular to which is a plane in which the material behaves in an 
elastically isotropic fashion. In the case of a unidirectional composite, the fiber axis is 
the preferred axis having a much larger elastic stiffness, and the isotropic plane is the 
plane perpendicular to the fiber axis with a much lower stiffness. 

For a transversely isotropic material there are five linear elastic (second order) 
moduli. This is the same as for the case of a hexagonal single crystal. However, there 
are only nine nonlinear (third order) moduli for transverse isotropy as opposed to either 
ten or twelve for the different point group symmetries of hexagonal single crystals. In 
the present work, the five linear (second order) elastic stiffness moduli for the unidirec- 
tional composites were computed from ultrasonic velocity measurements. These values 
were then compared with available data from the literature as well as checked for con- 
sistency with the assumption of transverse isotropy. From the linear elastic stiffness 
moduli, the compliance moduli were calculated and compared with those obtained by 
strain gauge measurements. These were also compared with reported values for a simi- 
lar material. Measurements were then made to determine the change in ’’natural” wave 
velocity as a function of applied stress and temperature. A variety of combinations of 
propagation direction and loading configurations were used with longitudinal and 
transverse waves. This data was also checked for consistency with the predictions of 
transverse isotropy. Then the velocity-stress data was used to calculate some of the 
nonlinear coefficients of this material. 
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H. LINEAR ELASTICITY 

H.A. Introduction 

In order to determine the nonlinear elastic moduli, the linear elastic stiffness and 
compliance moduli must first be measured. In this chapter, measurements of these pro- 
perties in a unidirectional graphite/epoxy composite will be discussed. Definitions of the 
terms used in linear elastic theory will be presented first, followed by a theoretical 
explanation of the methods used. Then, the measurements will be reported and com- 
pared with theoretically predicted and previously measured values for this material. 

II.B. Theory 

As discussed earlier, the theoretical basis for linear elasticity shows that it is 
only applicable for small or infinitesimal deformations. This permits the use of simplify- 
ing assumptions in the definition of strain for linear elasticity. The coordinates of a 
point of material in an undeformed body are defined to be a; with respect to the origin 
in an orthogonal coordinate system (Note: Einstein convention of summation over 
repeated indices is assumed throughout this paper). The same point of material defined 
by aj will be displaced to a new position Xj in the material after deformation. The dis- 
placement is a vector defined by 

Tf = Y - a* (2-1) 

having components Uj. The strain can now be defined in two ways. The Lagrangian 
strain tensor which provides a measure of the deformation with respect to the unde- 
formed material coordinates is defined by 

1 chi; <9llj <9u k <9u k 

_ 2 daj + <9a; + dajda. 


( 2 . 2 ) 
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In the Eulerian formulation of strain, the strain is described in reference to the 
deformed state and is given by 


e 


ij — 


1 

2 


5uj 5uj 

<9xj dxj 


5u k 5u k 

5a;5aj 


(2.3) 


The Lagrangian formulation is often referred to as a material description of strain while 
the Eulerian is called the spatial description. The Lagrangian is more often used in 
solid mechanics and the Eulerian in fluid mechanics. The assumption of infinitesimal 
deformation can now be applied. This assumption means that since 

5llj 5u; 

x - and x- «* (*•<> 


for infinitesimal deformations, then the quadratic term 


< 9 u k du k 


is much smaller 


{ d *\ d *\ J 

than the linear terms and can be dropped from the expression. Thus the small or 
infinitesimal strain tensor is defined by 


6ij 2 


duj duj 

daj ctaj 


(2.5) 


The small strain assumption also means that 

\ 


1 

C ij 2 


5u; chlj 


5 a 


5 a; 


1 _ 

2 


5uj 5uj 

5 x; 5 X; 


( 2 . 6 ) 


or that there is no difference between the Lagrangian and Eulerian formulations for 
strain. It can also be seen from the definition of the small strain tensor that it is a sym- 
metric tensor (i.e. ). 

The physical description of the components e n , and €33 is that they are the 
normal strains which provide the measure of the normalized change in length of the 
body along the respective coordinate axis. For example, if 1 0 is the unstrained length 
along axis Xj and lj is the strained length along X! , then 
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Oi - W 

e H “ i 


(2.7) 


The remaining strain components are the shear strains which are a measure of the 
change in angle between two originally orthogonal axes in the undeformed medium. 
This is demonstrated in Fig. 2.1 where the strain is given by 


€ 12 “ 


Axi 

x 2 


= tan (0) ~ 9 (for small strain). 


(2.8) 


The stress tensor (<7jj) provides a measure of the applied forces on the body. It is 
defined by 

^ _ force on the i ; th face in the j'th direction ^ 

1J cross sectional area of the i'th face 

A diagram of the nine components of the stress tensor is shown in Fig. 2.2. The 
engineering or nominal stress is used where the cross sectional area is that measured in 
the unstressed state. The true stress uses the instantaneous cross sectional area which 
is much more difficult to measure. The components cr n , and 033 are the normal 
stress components while the remainder are shear stress components. Laws of statics can 
be used to show that if a body is to be in rotational equilibrium, the stress tensor must 
also be symmetric. 

Having defined both stress and strain, the relationship between the two can now 
be discussed. For linear elasticity, the generalized Hooke’s law is used. This expresses 
the relationship between stress and strain as 

a \\ — c ijkl e kl ( 21 °) 

where Cy kl is the fourth rank linear elastic stiffness tensor. In the most general form this 
tensor has eighty-one components. However, because the stress and strain tensors are 
both symmetric, the following relations are valid 
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Figure 2.1 - Two dimensional representation of shear strain, 
solid line represents undeformed medium while 
dashed line is deformed medium. 
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c ijkl c jikl 


( 2 . 11 ) 


and 


c ijkl = c ijlk- 


( 2 . 12 ) 


This reduces the number of independent components to thirty-six. The assumed 
existence of a strain energy density function (<f>) which for the case of linear elasticity 
has the form 


= ^" c ijkl e ij e kl 


(2.13) 


provides further reduction of the number of independent elastic coefficients. It must be 
true that 

, 1 1 

g C ijkl e ij e kl g C ijkl e kl e ij- (2-14) 

Equating this expression term by term with respect to strain yields the relation 


c ijkl c klij- (2-15) 

This reduces the number of independent coefficients to twenty-one. At this point, the 
Voigt notation can be used to simplify things somewhat. For each pair of indices on 
stress, strain and elastic moduli, the following substitutions are made 

11 — *► 1 23 -»4 

22 -► 2 31 — ► 5 ( 2 . 16 ) 

33 — ► 3 12 — ► 6 


Using this notation, the generalized Hooke’s law can be rewritten as 


a A ~ c AB e B 


(2.17) 


where capital subscripts are summed from one to six. The elastic stiffness moduli can 
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now be written in matrix form for the most general anisotropic solid as 


C 11 

c 12 

c 13 

Cl4 

c 15 

c 16 

c 12 

c 22 

c 23 

c 24 

c 25 

c 26 

c 13 

c 23 

c 33 

c 34 

c 35 

c 36 

Cl4 

c 24 

c 34 

c 44 

C45 

c 46 

c 15 

c 25 

c 35 

C45 

C55 

c 56 

c 16 

c 26 

c 36 

c 46 

C 5 6 

c 66 


The reduction in the number of independent elastic constants can be carried 
further if material symmetries are taken into consideration. This is accomplished by 
examining the effect on the strain energy density function of a rotation of the material 
which places the material in an elastically equivalent condition due to symmetries of the 
material. Due to the scalar nature of the strain energy density function, it should have 
the same value in both configurations. An example of this procedure is now given. The 
strain energy density function is written out in terms of the twenty-one independent 
coefficients as 

<f) = ~(CnC? + C 22 c| + C 33 e 3 + C 44 C 4 + C 55 e 5 + c 66 e !) + c 12 e l e 2 

Jd 


+ c 13 e l e 3 + c 14 C l e 4 + c 15 e l e 5 + c 16 e l e 6 + c 23 e 2 e 3 + c 24 e 2 e 4 + c 25 e 2 e 5 
+ C 2 gC 2 Cg + C34C3C4 + C35C3C5 + C3gC3Cg + 0456465 + ^456465 + 0556565 ( 2 . 19 ) 


In an isotropic material any rotation is elastically equivalent. Assuming an isotropic 
material, a rotation of 180 degrees around the x 3 axis is chosen and is shown in Fig. 2 . 3 . 
The transformation matrix or direction cosine matrix [a^] is found by taking the cosine 
of the angle between the rotated and original axes. For this rotation it is given by 



-1 

0 

0 


0 0 

-1 0 

0 1 


( 2 . 20 ) 
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The effect of this rotation on the strain tensor is given by 

e ij = a ik aj,e k i. (221) 


Thus e ll — e H Cl Cl 

I 1 

^22 = e 22 ^2 e 2 

/ / 

e 33 = e 33 C 3 “ e 3 

, / ( 2 - 22 ) 

C 2 3 — e 23 or e 4 e 4 

I ■ 1 

e 31 = — e 31 ^5 e 5 

/ ' _ 

Cj2 — ^12 ^6 ^6 


The rotated strain energy density function can be written after substituting equation 
(2.22) as 

( f > = — (c u ef + C22C2 + 03363 4 - 04464 4 - 05565 4 - c 66 e|) 4 - 

2 

+ 0436^3 — 0446464 — 0456465 + 0456465 + 0336263 — 0246364 — 0356365 

+ C 2 6 e 2 e 6 ~~ c 34 e 3 e 4 — c 35 e 3 e 5 + c 36 c 3 e 6 + c 45 e 4 e 5 — c 46 e 4 e 6 — c 56 e 5 e 6)- ( 2 - 2 3) 

Evaluating term by term the relation 

0' = (j) (2-24) 


yields 


c 14 — c 14 

C 15 = ~ c 15 
c 24 = “ c 24 
025 = — c 25 
O34 = *“ c 34 
c 35 — “ c 35 
C46 = — c 46 
c 56 = — c 56 
all other Cy = Cy. 


(2.25) 


If equations (2.25) are to be true, then 
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c 14 ~ c 15 “ c 24 — c 25 — c 34 = c 35 ~ c 46 = c 56 = ( 2 - 26 ) 

Application of further symmetry operations will show that for an isotropic material, the 
number of independent elastic coefficients reduces to two. The matrix then appears as 




(2.27) 


For a transversely isotropic material, the number of independent coefficients is five while 
for an orthotropic material it is nine. The derivation of the independent coefficients as 
well as their final matrix forms for these two elastic symmetries which are important for 
unidirectional fiber reinforced composites is shown in Appendix A. 


The relationship between strain and stress is given by 

^ij = s ijkl <:r kl * (2.28) 

The moduli form the linear elastic compliance tensor. The conditions of symmetry 
for the stiffness tensor also apply to the compliance tensor. The relations between the 
compliance and stiffness moduli for a given symmetry can be determined by inverting 
the stiffness matrix. 

The value of the elastic moduli are also dependent upon the conditions (isother- 
mal or adiabatic) in which they are measured. Isothermal moduli are measured at con- 
stant temperature (T) while adiabatic moduli are measured at constant entropy (S). 
The work by Brugger [12] provides a precise thermodynamic definition of the elastic 
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moduli and the relationship between adiabatic and isothermal coefficients. In terms of 
the free energy (F), the isothermal stiffness moduli are defined by 


T 

c ijkl = Po 


& F 


dViidvu 


(2.29) 


where p Q is the unstrained density. 
The adiabatic stiffness moduli are 


S 

c ijkl = Po 


&U 


{dviidvu} 


(2.30) 


where U is the internal energy. 

The adiabatic and isothermal compliances are defined in terms of the enthalpy (H) and 
the Gibbs function (G) respectively and the thermodynamic tensions (t^) which are 


t’ij ~ PO 


f dU ' 


f 

d F 

#1 

= Po 

5 



(2.31) 


The compliance moduli are 


s 







and 







(2.32) 


(2.33) 


The relationship between adiabatic and isothermal compliance moduli was given in 
matrix form by Hankey and Schule [40] as 


T S 
S AB = S AB + 


[P 0 CpJ 


f*A a B 


(2.34) 


where c p is the specific heat at constant pressure, 
and a A is the linear coefficient of thermal expansion defined by 
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a A = 


d*/A 
dT ■ 


(2.35) 


The difference between adiabatic and isothermal moduli is usually quite small and often 
less than experimental uncertainty. 

The determination of the adiabatic linear elastic stiffness moduli can be made 
by measuring the density and the ultrasonic wave velocity. The derivation of the equa- 
tions of motion for the propagation of linear elastic waves in an anisotropic medium 
provides the theoretical basis for this technique. This subject has been treated by a 
number of authors among the earliest being Love [41] and it is reviewed in detail by 
Green [5]. In deriving the wave equations, a number of initial assumptions about the 
propagation medium are made. It is assumed to be unbounded, homogeneous and con- 
tinuous. Next, an infinitesimal volume element of the medium is examined and the net 
unbalanced forces acting on it are determined. This is accomplished by finding the 
resultant of the variation of stress across the element for each direction and multiplying 
this by the respective cross sectional area. Fig. 2.4 shows the variation of stress along 
the Xi direction for the infinitesimal volume element with sides of length Sx l} foc 2 and fa 3 . 
If there were no variation of stress and thus no net unbalanced forces, there would be 
no wave motion. The net unbalanced force along x { is given by 


[(°n + 


dcr n 

-q— 1 &i) - <r n ]6x 2 8x 3 + [(<r 21 + 


dcr 2l 

- q — 6 * 2 ) ~ 02i]&i&3 


+ [(^l + 


dao 1 

(5x 3 ) - <7 31 ]&C,&C2. 

3 


(2.36) 


This can be reduced to 


da. 


11 


<9xi 


+ 


da 


21 


dx 0 


+ 


da 


31 


dx.o 


|<5x 1 <5x 2 • 


(2.37) 


Similar expression can be written for the net unbalanced force along x 2 and x 3 . Accord- 
ing to Newton’s second law these unbalanced forces must equal rniij where m is the mass 
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of the element and \i { are the second time derivatives of the displacements which are the 
acceleration components. Actually, to be more exact, it is pointed out that Newton’s 
second law states that the force is equal to the time derivative of the momentum pj . 
However, 

Pi ~ mUj (2.38) 


and therefore 


Pi = mUj + muj. 


(2.39) 


Since the mass of the element is assumed to be constant in time the forces are simply 
given by muj . Therefore the balance of forces can be written in component form as 


mUj = 


dcr V] 


dx 


J ) 


kx 1 &C 2 & 3- 


(2.40) 


It is noted that the body forces are neglected in this equation. The mass can be rewrit- 
ten in terms of the density p 0 and the volume (V) which is given by 

V =* &C!&C 2 &C 3 (2.41) 


m — PqSx^Sx. 2<5xg. 


(2.42) 


Substituting (2.42) into (2.40) and canceling appropriate terms yields 

Po“i = 




dX; 


I ) 


(2.43) 


To rewrite this in terms of only density and displacement, the stress is first rewritten in 
terms of the strain using Hooke’s law (2.10). Then the strain definition (2.6) is used to 
produce the expression for stress in terms of displacement as 




Now that the equation of motion has been derived, the next step is to choose a 
solution for the displacements as a function of time and space. For a plane wave, a 
solution can be written as 

( 2 - 48 ) 

where Ao is the amplitude of the wave, 

are the direction cosines of the particle displacement vector, 
t is the time, 

c j is the angular frequency 
and k m are the components of the wave vector. 


The wave vector components (k m ) are related to the direction cosines of the wave nor- 
mal (l m ) and the wavelength (X) by the expression 
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v / 


(2.49) 


where k is the wave number. 

The solution for the displacement can now be substituted into the equation of motion. 
This will lead to the following expression 

Cijklk 2 lllj<* k = Po<J<*i (2.50) 

which can be rewritten as 

(CijklVj - A)V 2 %K = 0 (2.51) 


where v is the phase velocity of the elastic wave and is given by 


V = 


u i 

k ‘ 


(2.52) 


In equation (2.51), it is noted that o k is arbitrary and therefore not necessarily equal to 
zero. In order for this equation to have nontrivial solutions, it must be true that the 
determinant of the matrix of coefficients must equal zero. That is 

IcijklMj -Po v2 $k| =0- (2.53) 

To more conveniently write this equation, a matrix [X ik ] is defined by 

\k 2=8 c ijklMj- (2-54) 

Therefore, equation (2.53) becomes 

l X ik -Po y2<5 ik| =0. (2.55) 

This is a characteristic equation which when expanded forms a cubic expression in terms 
of p 0 v 2 . Since X ik is a symmetric matrix, the three solutions for p 0 v 2 must all be real. 
The physical significance of this is that for any arbitrary direction in an anisotropic 
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medium, three plane waves may be propagated. The velocity of each of these waves 
depends on the density, the elastic constants and the direction cosines of the wave nor- 
mal. Thus, if the ultrasonic wave velocities are measured for a number of directions in 
an anisotropic medium with known density and elastic symmetry, the linear elastic 
stiffness moduli can be determined. Although not formally designated as such in the 
preceding equations, these moduli will be the adiabatic stiffnesses. 

The particle displacement vector direction cosines can be determined for each of 
the three waves propagating along a given direction. This is accomplished by solving 
for the eigenvectors for each eigenvalue (p 0 v 2 ) using the equation 

pQ\ 2 a x - \ ik a k , (2.56) 

and also applying the fact that 

a x a { = 1. (2.57) 

These direction cosines determine the mode of the wave. If they are the same as the 
wave normal direction cosines, that is 


<*i = lj, (2.58) 

then the particle displacement of the elastic wave lies along the direction of propaga- 
tion of the wave. A wave of this type is referred to as a pure mode longitudinal 
(compressional) wave. If the particle displacement is perpendicular to the propagation 
direction and therefore 


Ojli = o, (2.59) 

then the wave is called a pure mode transverse (shear) wave. Illustrations of longitudi- 
nal and transverse mode waves are given in Figs. 2.5 a) and b). 
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Figure 2.5 - (a) Illustration of a longitudinal (compressional) 

elastic wave demonstrating directions of particle 
vibration and propagation. 
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Figure 2.5 - (b) Illustration of a transverse (shear) wave 



- 38 - 


However, in general the wave need not be pure mode. It can be either quasilong- 
itudinal or quasitransverse depending on which type of particle motion is more dom- 
inant. Of the three waves propagating along a given direction, one is longitudinal or 
quasilongitudinal while the other two are transverse or quasitransverse. 

Another consideration of linear elastic wave propagation in an anisotropic 
medium is that the direction of the energy flux of the wave may not in general be the 
same as the direction of the wave front normal. That is, there may be a refraction of 
the wave because of the anisotropy of the medium. The energy flux vector is defined as 
the direction of the flow of energy per unit time per unit area. A diagram of the energy 
flux deviation of a quasitransverse wave is shown in Fig. 2.6. Equations for expressing 
the energy flux vector and its deviation from the direction of the wave normal are 
reviewed by Green [5] and are presented here. 

First, the expressions for the kinetic energy (K) and the potential energy (P) 
which is contained in the wave field are given by 

K = —J puj 2 dr (2.60) 

£ 

and 

p = jJV ij dr ( 2 - 61 > 

where dr is a volume element of the wave field. 

The sum of K and P yields the total elastic energy E contained in the wave field which 

is 

= = j (/puj 2 dr + JVijCijdrj. (2.62) 

The time derivative of E is then given by 


(2.63) 
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Figure 2.6 - Diagram of a quasitransverse elastic wave undergoing 
energy flux deviation illustrating the directions of 
particle displacements, wave normal, and energy flux 
vector. 
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To simplify this expression the term under the second integral is now considered. It can 
be rewritten as 



cr ij e ij 


+ 


(2.64) 


d , v . 

Hooke’s law (2.10) can now be substituted into this equation to express ^ in 

terms of the elastic moduli and the strains as 


(2 

■^( <T ij € ij) " c ijkl*k! e ij + c ijkl e kl e ij- 


(2.65) 


Because of the symmetry of the elastic moduli, equation (2.65) can be rewritten as 



c ij) — c klij^kl e ij c ijkl e kl e ij- 


( 2 . 66 ) 


Again using Hooke’s law this becomes 


(°ij e ij) — ^kl^kl + 


(2.67) 


which can be rewritten as 




( 2 . 68 ) 


after the dummy indices k and 1 are summed out of the equation. The time derivative 
of the strain in equation (2.68) can be rewritten in terms of the displacements using the 


definition of strain as 


1 duj diij 

2 dx\ dxj 


(2.69) 


This is substituted into equation (2.68) resulting in 



(9Uj dUj 

<9Xj <9x; 


(2.70) 
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which because of the symmetry of the stress tensor is 

d <9u; 


(2-71) 


Equation (2.71) is now put back into equation (2.63) to give the time derivative of the 
total energy of the elastic wave field as 


h (^ 11 : 

ir = M u i dr + M^r dr - 


(2.72) 


a , . . , . wji 

dx-, (l7 ‘ iU|) aii dx> + u ' ax. ’ 


(2.73) 


d , .v 

a » X" = •3TK u i) ~ «iX" 


dxj <9xj 


(2.74) 


the second integral in equation (2.72) can be rewritten as 

f i r_3 , • x, r. , 

/ < 7 ii' 3 r dr “ /^r( CT u u i) dT - / u ,^r d ’'- 


(2-75) 


Gauss’s theorem can now be applied to the first integral of equation (2.75) to rewrite it 
as a surface integral of the form 


Jxrfa i^) dr = 


(2.76) 


This can be substituted into equation (2.75) which is then inserted into equation (2.72) 
which gives 


rP 3(7; ; I 

aT - / ^-^-M 7 - + Kii d s. 


(2.77) 
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The first integral equals zero because of the equation of motion (2.47). Therefore, it 
remains that 

If- - Ku(dS. (278) 


This equation represents the energy flow into the volume element r through the surface 
S. The energy flux vector with components Ej can be formed which gives the energy flux 
flowing outward from the volume r through S by 

Ej = - CTjjUj. (2.79) 

If the time averaged values E of the energy flux vector over one cycle of the elastic 
wave are computed and it is defined that 

E| = Ei + E| + E|, (2.80) 


then the direction cosines of the energy flux vector (E[) can be formed by 



(2.81) 


If 



(2.82) 


then the wave propagates along the direction of the wave normal and there is no energy 
flux deviation. Otherwise, the angle ( 6 ) of deviation of the energy flux from the wave 
normal can be determined by 

S = COS -1 (Ejl i ). (2.83) 


Therefore, to determine the elastic moduli by using ultrasonic velocity measure- 
ments, a number of considerations must be examined. First, the elastic symmetry and 
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thus the form of the elastic stiffness matrix must be known. Then a direction of propa- 
gation must be chosen and the solutions for the velocities of the three waves that pro- 
pagate along this direction must be calculated in terms of the density and the elastic 
moduli. Also, any energy flux deviation of the three wave modes must be determined. 
This is repeated for a number of different directions until enough independent equations 
are derived to calculate all of the elastic moduli from the corresponding velocity meas- 
urements. Usually, an attempt is usually made to choose directions which propagate 
waves which are pure mode and suffer no energy flux deviations to make the measure- 
ments as easy as possible. Unfortunately, this can not always be accomplished. 
Brugger [42] published an equation which can be used to determine the pure mode direc- 
tions for the various elastic symmetries. In general form this equation is 

e ijk c jlrs l kl|Ms = 0 ( 2 - 84 ) 


where € V)k is the permutation tensor. 


It is defined by 



if any i,j,k are equal 
if i,j,k are in cyclic order 
if i,j,k are in acyclic order 


(2.85) 


Also presented are the pure mode wave directions for the various elastic symmetries. 

The specific case of elastic wave propagation in a unidirectionally graphite fiber 
reinforced composite will now be discussed. First, the applicability of the previously 
derived wave equations to a composite material must be considered with respect to the 
heterogeneous nature of the medium. This violates one of the initial assumptions that 
the medium be homogeneous. However, it is generally assumed that if the wavelength 
of the elastic waves is much larger than the size of the inhomogeneity, then the material 
is macroscopically homogeneous and the wave equations are applicable. In the materi- 
als studied, the size of the inhomogeneity is given by the diameter of the fibers which is 
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on the order of seven //m . In the frequency range used (« 2.25MHz) , the wavelengths 
for the different wave modes and directions of propagation ranged from about 700/im to 
almost 5000/im . Therefore, the previously stated assumptions are valid as are the wave 
equations. 

Another important question is the use of Hooke’s law because of the viscoelastic 
nature of most polymeric materials. However, epoxies and their composites, because of 
the high degree of crosslinking along the molecular chains, generally exhibit little viscoe- 
lastic behavior. This is especially true at the temperatures at which these measure- 
ments were made (room temperature). Measurements of the wave speeds and thus the 
moduli as a function of frequency demonstrated no viscoelastic effects. Therefore the 
application of Hooke’s law appeared to be valid. 

The appropriate elastic symmetry model for a unidirectionally fiber reinforced 
graphite/epoxy composite must now be chosen in order to derive the specific solutions 
for the velocities of ultrasonic wave propagation. The most common symmetry model 
used for unidirectional composities is that of transverse isotropy. In a transversely iso- 
tropic material, there exists a preferred axis perpendicular to which there exists a plane 
in which the material behaves isotropically. In unidirectional graphite/epoxy, the fiber 
direction (x 3 ) is the preferred axis while the x,x 2 plane is the plane of isotropy as shown 
in Fig. 2.7. Transversely isotropic behavior of unidirectional composites is based on the 
assumption of a random distribution of fibers in the XjX 3 plane as well as alignment of 
the fibers along the x 3 axis. The symmetry conditions for such a material are that it be 
elastically insensitive to twofold (180 degree) rotations about each of the three axes as 
well as insensitive to any rotation around the x 3 axis. Application of these symmetry 
conditions as shown in Appendix A demonstrates that there are only five independent 
elastic moduli. The matrix form of the elastic stiffness moduli for a transversely isotro- 
pic material is then given by 



Illustration of axis designation in a unidirectional 
fiber reinforced composite with respect to the fiber 
direction and the lamina stacking direction. 
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[ c ab1 — 


C 11 

c 12 

c 13 

0 

0 

0 

c 12 

Cll 

c 13 

0 

0 

0 

c 13 

c 13 

c 33 

0 

0 

0 

0 

0 

0 

C44 

0 

0 

0 

0 

0 

0 

c 44 

0 

0 

0 

0 

0 

0 

Cll - 

2 


( 2 - 86 ) 


However, if for some reason the fibers are not randomly spaced in the X 1 X 2 plane, 
this plane will no longer be elastically isotropic. In this case, the material will possess 
an orthotropic symmetry. This might occur if because of the lay up process of manufac- 
ture there was a resin rich layer between each lamina as shown in Fig. 2.8 or if some of 
the laminae were misoriented. The resulting orthotropic material would now have only 
the two fold symmetries about each of the three axes as is the case with cross ply or 
angular ply composite lay ups. For this symmetry the number of independent linear 
elastic moduli increases to nine as shown in Appendix A. The matrix form of the moduli 
is then 


[ c ab1 = 


C 11 

c 12 

c 13 

0 

0 

0 


c 12 

c 22 

c 23 

0 

0 

0 


c 13 

c 23 

c 33 

0 

0 

0 


0 

0 

0 

c 44 

0 

0 


0 

0 

0 

0 

C55 

0 


0 

0 

0 

0 

0 

c 66 


( 2 . 87 ) 


Ultrasonic velocity data can be compared with the equations derived for both 
models to determine which gives the best fit to the data. The derivations of the equa- 
tions for the velocity as a function of density and elastic moduli for a number of 
different directions and pure wave modes for transversely isotropic and orthotropic sym- 
metries are presented in Appendix B. However a summary of the results is presented in 
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Figure 2.8 - Illustration of unidirectional composite with resin 
rich regions between laminae causing the material to 
have orthotropic elastic symmetry. 
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Table 2.1. To determine all of the elastic moduli for the two models, several nonpure 
mode waves must be used. The results for these waves are presented in Table 2.2. 


II.C. Experiment 

In most of the ultrasonic velocity measurements made in this study, a pulse echo 
overlap velocity measurement system was used. This system was briefly discussed in the 
Introduction. However, the details of its operation will now be discussed. Figure 2.9 
shows a block diagram of the entire pulse echo overlap system. A Matec model 6600 
pulse modulator and receiver was used to generate an electrical tone burst of variable 
amplitude and frequency. This was used to excite an ultrasonic transducer (typically 
undamped crystals of quartz or PZT with fundamental resonances of approximately 
2.25 MHz.) which excited ultrasonic waves into the graphite/epoxy samples. The result- 
ing echoes from the back surface of the specimen were detected by the same transducer 
and converted back into electrical signals. The received signals were input into the 
receiver and amplified. After amplification, they were displayed on a Hewlitt Packard 
(H.P.) 1743A dual trace oscilloscope. A typical pulse echo pattern is shown in Fig. 2.10. 
The echoes appear identical to the generating signal (main bang) except for an exponen- 
tial decrease in amplitude due to attenuation in the sample. 

To have the pulse echo pattern appear stable on the oscilloscope, the Matec gen- 
erator and the H.P. 1743A oscilloscope must be triggered simultaneously. This was 
accomplished by using the triggering components of the pulse echo overlap system. The 
signal from a continuous wave (c.w.) oscillator (Matec model 110) was input into a 
Matec 122B decade divider and dual delay generator. The Matec 122B divided the ori- 
ginal signal by a preset amount (typically 1000) and then generated a trigger signal at 
the divided frequency. This is demonstrated for frequency division by five in Fig. 2.11. 
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Table 2.1 - Equations for pure mode linear elastic wave propagation in 
transversely isotropic and orthotropic media. 
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Table 2.2 - Equations for non pure mode, off axis linear elastic wave 
propagation in transversely istotropic and orthotropic media. 


Transverse Isotropy 


Propagation Wave Po y2 

Direction Mode 

Cosines 



lj = 0 QL Same as QL above 

1 2 = PT Same as PT above 


1 3 = QT Same as QT above 



PL 



la = 0 


PT 


j(cn - C12) 
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Figure 2.9 - Block diagram of pulse echo overlap system. 
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Figure 2,10 - Illustration of a typical pulse echo pattern 
demonstrating voltage versus time. 
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The trigger signal was then used to trigger both the Matec model 6600 and the oscillo- 
scope. The Matec model 122B also generated a dual strobe which was input into the z 
axis of the scope. This had the effect of intensifying two chosen regions of the oscillo- 
scope trace. The width and spacing of the strobes was variable. 

The high accuracy in measuring the time between two echoes needed to compute 
the velocity arises because of a special triggering technique. First, the echo pattern was 
set up as previously described. The strobes were set so that two echoes were intensified 
on the trace and the scope intensity control was adjusted so that only these two echoes 
remained visible on the trace. Then, using another output of the Matec model 122B, a 
trigger signal was applied to the oscilloscope which was at the original rate of the 
Matec oscillator which was typically 1000 times faster than the trigger rate of the 
Matec model 6600. This caused an oscilloscope trace to be drawn out 1000 times during 
the single pulse echo cycle of the Matec 6600. The effect of this, because of the per- 
sistence of the oscilloscope screen, was to make the two echoes appear together on the 
screen as shown in Fig. 2.12. If the trigger frequency was adjusted, the two echoes could 
be caused to coincide on the screen or "overlapped" as shown in Fig. 2.13. The fre- 
quency (f) of the Matec c.w. oscillator at which this occurred was related to the time 
between echoes (t) by 

t = j (2.88) 


where n was an integer. 

By measuring two successive overlap frequencies the integer n could be taken out 

of the equation and the time of flight then given by 


t = 



(2.89) 


i 


A number of overlap frequencies were measured and the results averaged to get a more 
reliable value for t. The main uncertainties of the measurement of t were the 



-55- 


C. W. Source 



A A. A JL 

Output of Motec 122 8 


Figure 2.11 - Illustration of frequency division of a continuous 
wave (c.w.) source by the Matec model 122B by a 
factor of five. 


Figure 2.12 - Illustration of an oscilloscope trace of two echoes 
during the pulse echo overlap technique where the 
echoes are not properly overlapped. 



Figure 2.13 - Illustration of two properly overlapped echoes. 
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uncertainty in the operator’s visual ability to detect the overlapped condition and the 
uncertainty in the measurement of the c.w. oscillator frequency. 

Because of time delays of the ultrasonic pulse in the bond and phase shifts at the 
specimen-transducer interface, corrections had to be made to determine the true time of 
flight (0) from the measured time (t). The theory for these corrections was presented by 
McSkimmin [19] and is outlined in Appendix C. In this theory it is necessary to deter- 
mine the frequency of the ultrasonic pulse. This was accomplished by inputting the tone 
burst into a mixer (Anzac MD-40) along with the signal from a HP. 3325A frequency 
synthesizer/function generator. The output of the mixer was then input into a Tek- 
tronix 53C oscilloscope. The output of the mixer was the sum and difference frequencies 
of the tone burst and the HP. c.w. source. If they were of the same frequency, the out- 
put of the mixer was a constant voltage (0 Hz) dependent on the phase difference 
between the two signals. This phase difference was variable each time the Matec 6600 
was triggered and thus a series of parallel lines was generated on the scope trace. This 
is shown in Fig. 2.14. If the frequencies were different, then the output would have the 
difference frequency again phase shifted on each new trace. Therefore, the frequency of 
the tone burst from the Matec 6600 could be monitored by the frequency of the HP. fre- 
quency synthesizer. 

The remaining equipment was used to collect data. A Tektronix 4051 computer 
was used to read the trigger frequency values and calculate the time of flight which was 
used to calculate the velocities and the moduli. The program used to control the pulse 
echo overlap data acquisition and analysis is shown in Appendix D. 

There were a number of factors which affected the size and geometry of the sam- 
ples used. They had to have parallel sides to be used in the pulse echo technique. Also, 
the thicker the specimen, the more accurately its time of flight could be measured. 
However, because of the high ultrasonic attenuation of graphite/epoxy composites, espe- 
cially for shear waves, there were limits on the maximum thickness. Also, limits on 
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Figure 2.14 - Illustration of the oscilloscope trace of the output 
of a mixer when the input c.w. and pulse sources are 
of the same frequency. 


| 
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sample thickness occurred because of the desire to have the wavefront remain planar. 
The distance (d) over which the wavefront remains planar was given by Mason [43] as 



where R is the radius of the transducer 
and X is the wavelength of the ultrasonic wave. 

However limitations on the thickness because of the large attenuation usually occurred 
before this requirement became important. There were also restrictions on the lateral 
dimensions of the samples. The equation of motion are based on the assumption of a 
wave propagating through a unbounded medium. However, ultrasonic velocities have 
been shown by Tu et al. [44] to be independent of lateral dimensions if the sample was 
at least two and one half times the acoustic wavelength. 

The above constraints provided the basis for the choice of the sample geometries 
used. Ideally, the measurements of the elastic moduli would be made using only one 
sample to avoid any problems with sample to sample variations in the elastic properties. 
These problems are especially prevalent in polymer composite materials where elastic 
property variations may occur not only from batch to batch due to different curing pro- 
cedures but also within a given batch because of regional variation of fiber content, void 
content and other important parameters. However, because of the above constraints 
and the necessity to use off axis non pure mode wave propagation to determine some of 
the moduli, several samples were required. These samples were all cut from a 150 ply 
laminate of T300/5208 unidirectional graphite/epoxy that was nominally 12 X 15 X 0.9 
in. (30.5 X 38.1 X 2.3 cm.) after cure. This panel was designated MCIS 1 and layed up 
and cured at NASA, Langley Research Center. The finished panel was ultrasonically 
C-scanned by the Materials Characterization and Instrumentation Section at NASA 
Langley for the presence of gross voids or delaminations. None were found to be 
present. A number of tensile and compression samples were then cut and machined 
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from this panel. Those used for the ultrasonic testing were MCIS 1.10A, 1.8A, 1.10B, 
1.7A, 1.5 A,B>E and F and 1.6A. whose dimensions are listed in Tables 2.3 and 2.4. The 
samples were generally flat and parallel to +/- 0.0003 in. (0.00076 cm.) or better. Sam- 
ples MCIS 1.10A, 1.8A, 1.10B and 1.7A were all rectangular parallelepipeds with their 
faces oriented perpendicular to the three coordinate axes. The densities for each are 
also given in Table 2.3 which show some sample to sample variation. The remaining 
samples were cut to allow wave propagation in a direction that was 45 degrees between 
two of the axes and perpendicular to the third. These samples were used to propagate 
the waves listed in Table 2.2. The sample dimensions are given in Table 2.4. The densi- 
ties of these samples were not measured but the value of 1.54 g/cm 3 was used for all cal- 
culations involving these specimen. 

Measurements of the pure mode longitudinal wave speeds along each of the three 
coordinate axes were made using specimen MCIS 1.10A. Specimen MCIS 1.8A, 1.10B 
and 1.7A were used to determine the six different combinations of polarization and pro- 
pagation directions for the pure mode shear waves along the three axes. The velocities 
and the product p 0 v 2 for each of these waves are shown in Table 2.5. The uncertainties 
for the shear wave velocities were larger than those of the longitudinal velocities 
because of the much higher attenuation of shear waves. Even using the thinner speci- 
men for the shear waves did not overcome the high attenuation. This resulted in a 
much smaller amplitude second echo which was difficult to reproducibly overlap. The 
uncertainty in the longitudinal velocity along x 3 was also large because of the high velo- 
city and thus the short time of flight along this direction. 

The measured values were then compared with theory to check for transversely 
isotropic behavior. The products p 0 v 2 calculated using the longitudinal waves along 
Xi and x 2 in sample MCIS 1.10A should be equal as should the products calculated from 
the two shear waves propagating along x 3 in sample MCIS 1.8A. Comparisons for the 
other similar products were not valid because they were obtained on different samples. 
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Table 2.3 - Sample dimensions and densities for samples used for pure 
mode elastic wave propagation. 


Sample 

Number 

Dimensions (in.) 

Xj x 2 x 3 

Density 

(g/cm 3 ) 

1.10A 

0.8001 

0.7999 

0.8002 

1.5404 +/- 0.0004 

1.10B 




1.5460 +/- 0.0006 

1.8A 

0.8002 

0.8000 


1.5384 +/- 0.0005 

1.7A 

0.5001 

0.8001 

0.7999 

1.5411 +/- 0.0005 


Table 2.4 - Sample dimensions for samples used for off axis, non pure 
mode wave propagation. 


Sample 

Number 

1.5A 

1.5B 

1-5E 

1.5F 

1.6A 


Dimension (in.) along axis 
45 degrees between axes 


Dimension (in.) along axis 


0.4499 

0.5000 

0.4491 

0.4993 

0.5001 


*i, x 3 

Xl> x 3 

*2, X 3 

X 2 , X 3 

Xi, x 2 


0.9 

0.9 

0.9 

0.9 

0.9 


Xg 

X2 

h. 

X 3 


Table 2.5 - Data for pure mode linear elastic wave propagation. 


Sample 

Number 

1.10A 

1-10A 

1.10A 

1.7A 

1.7 A 

1.10B 

1.10B 

1.8A 

1.8A 


Propagation Polarization 

Direction Direction 


Xl 


*1 


x 2 


x 2 


X 3 


X 3 


X 2 


h. 

Xg 

X2 

X3 

X 3 


Xs 

ill 

^3 

Xl_ 

X 2 


Velocity 

(l0 5 cm/s) 

3.0467 +/- 0,0005 

3.0390 +/- 0,0008 

8.390 +/- 0,007 

1.557 +/- 0,002 

1.849 +/- 0,003 

1.558 +/- 0.003 

1.850 +/- 0,003 

1.850 +/- 0,002 

1.851 +/- 0.002 


/> 0 v 2 

(GPa) 

14.295 +/- 0.005 
14.226 +/- 0.006 
108.4 +/- 0.1 
3.736 +/- 0,007 
5.27 +/- 0,01 
3.75 +/- 0,01 
5.29 +/- 0,01 
5.265 +/- 0.008 
5.271 +/- 0.008 

















































- 62 - 


In the first comparison, the difference was larger than experimental error. However, it 
was still quite small and thus the deviation from transversely isotropic elastic behavior 
was small. In the comparison of the shear wave products, the difference is smaller than 
experimental uncertainty. Thus, overall the elastic behavior of the material was very 
nearly transversely isotropic for linear elasticity. 

Using the data presented in Table 2.5, four of the transversely isotropic moduli 
can be calculated using the equations from Table 2.1. Also, six of the nine orthotropic 
moduli can be calculated. However, to determine the transversely isotropic c 13 and the 
orthotropic c 12 ,Ci 3 , and c^ , off axis quasilongitudinal and quasishear mode waves had to 
be used. The samples used for these measurements were MCIS 1.5 A,B,E and F and 
1.6A. As pointed out by Kriz and Stinchcomb [45], these wave modes suffer energy flux 
deviations as large as 43 degrees which impose further specimen size requirements. The 
basic requirements are that the sample was thick enough to allow the two quasi-mode 
waves to separate into distinct wave forms and thin enough so that the deviation of the 
wave did not cause it to impinge on the side of the specimen. The energy flux deviation 
and large attenuation also prevented the use of the pulse echo method for the velocity 
measurements. For these measurements, two transducers were used in a through 
transmission arrangement. The time of flight was first measured through a fused quartz 
delay line. Then it was measured through the delay line and the specimen for the 
different wave modes. The difference between the two was taken to be the time of flight 
and used to calculate the velocity. The uncertainty of this method is much larger than 
the pulse echo overlap method resulting in much reduced accuracy. The results are 
shown in Table 2.6 which show large scatter due to uncertainty in the measurement as 
well as material variations. 

Using the data from Tables 2.5 and 2.6 with the equations from Tables 2.1 and 
2.2, the independent elastic moduli were determined for both the transverse isotropy 
and orthotropy models. The results are presented in Table 2.7. Where possible, the 
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multiple results were averaged except for c 13 and C23 where there was a wide range of 
values. This wide range was due to the large uncertainty in the measurements which 
propagates into an even larger uncertainty in the calculated values. The values 
presented for c 13 and c^ were calculated using the quasishear data because the quasi- 
longitudinal data yielded complex values for the moduli. 


The values for the transversely isotropic model can be compared with those 
reported by Kriz and Stinchcomb [45]. They reported data for the moduli of T300/5208 
graphite/epoxy measured ultrasonically as well as the values theoretically predicted 
from the properties of the fibers and matrix. The theoretical values are derived for an 
assumed volume fraction of fibers of 0.67. These are all presented in Table 2.8. The 
values all compare favorably except for c^ which was much lower in the present study. 
The lower density of the samples used in the present experiments implies that there may 
be a lower fiber content and/or a higher void content. This would account for some of 
the discrepancy in the value of c^ . There may also be some misalignment of the fibers 
along the x 3 axis which would also lower the value of c^ although there is no other 
confirming evidence of this fact. The values of Kriz and Stinchcomb also show much 
more scatter demonstrating sample to sample variations. 


To calculate the nonlinear elastic moduli, the isothermal linear elastic compli- 
ance moduli are also needed. The adiabatic compliance moduli can be obtained by 
inverting the stiffness matrix. Corrections for differences in adiabatic versus isothermal 
moduli can then be applied if they are larger than the uncertainties of the inverted 
values. For the case of transverse isotropy the compliances are given by the following 
equations, 


S 11 “ 


c ll c 33 — c 13 


c 33( c ll C 12) + 2c 1 2 3 (c 12 — C n ) 


(2.91) 


C 11 + c 12 


s 33 “ 


c 33( c ll + c 12) 2cf 3 


(2.92) 
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Table 2.0 - Data for off axis, non pure mode linear elastic wave propa- 
gation. 


Sample 

Number 

Prop 

Diri 

agation 

action 

Wave 

Mode 

Velocity 

(10 6 cm/s) 

P ov 2 
(GPa) 

Assuming 
p 0 = 1.54g/cm 3 

Perpendicular to 

45 degrees between 

1.5A 

x 2 

Xi, x 3 

QL 

5.78 

51.4 

1.5A 

x 2 

Xi, x 3 

QT 

2.47 

9.40 

1.5B 

X 2 

*1, X 3 

QL 

5.73 

50.6 

1.5B 

x 2 

Xj, x 3 

QT 

2.38 

8.72 

1.5E 

Xl 

x 2 , X 3 

QL 

6.00 

55.44 

1.5E 

*1 

x 2 , x 3 

QT 

2.42 

9.01 

1.5F 

*1 

x 2 , X 3 

QL 

5.86 

52.9 

1.5F 

X 1 

x 2 , X 3 

QT 

2.40 

8.87 

1.6A 

X 3 

Xi, X 2 

QL 

3.043 

14.26 

1.6A 

*3 

Xi, X 2 

QT 

1.558 

3.74 


Table 2.7 - Measured linear elastic moduli for transversely isotropic 
and orthotropic symmetries. 


Transversely 

Value 

Orthotropic 

Value 

Isotropic 

Modulus 

(GPa) 

Modulus 

(GPa) 

cn 

14.26 

Cll 

14.295 

c 12 

6.78 

c 12 

6.78 

c 13 

3.0 - 8.9 

c 13 

3.3 - 9.0 

c 33 

108.4 

c 33 

108.4 
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Table 2.8 - Comparison of measured transversely isotropic stiffness 
data with that presented by Kriz and Stinchcomb [45]. 


Modulus 

Present 

Work 

(GPa) 

Kriz and Stinchcomb 

i 


Experiment 

Theory 

c n 

14.26 

15.0- 15.7 

14.5 

c 12 

6.78 

7.10 - 7.58 

7.24 

c 13 

3.0 - 8.9 

6.96 - 9.09 

6.50 

c 33 

108.4 

154 

161 


5.27 

7.84 

7.10 
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— ( c 12 c 33 — cf 3 ) 

(2.93) 

C 33( C U ~ c 1 2 2) + 2c 1 2 3( c 12 — c ll) 

c 13( c 12 — c ll) 

(2.94) 

c 33( c n — c n) + 2 c 1 2 3 (c 12 — c n ) 


and 



(2.95) 


Because of the large uncertainty in c l3 , there is a variation in the possible values for 
the compliance moduli. The values for the compliances for the two extreme values of 
Cjg are given in Table 2.9. The wide variation in the values of c 13 has little effect except 
on the modulus s 13 . 

The inversion of the orthotropic stiffness moduli yields the following equations 


c 22 c 33 — c 23 


S 11 — 


c ll( c 22 c 33 — c 23) — c 12( c 12 c 33 — 2c 13 c 23) c 13 c 22 


(2.96) 


s 22 — 


s 33 _ 


c ll c 33 — c 13 


c ll( c 22 c 33 — c 23) — c 12( c 12 c 33 — 2c 13 c 23) C 13 C 22 

2 

c ll c 22 ~ c 12 

c ll( c 22 c 33 — c 23) ~ c 12( c 12 c 33 — 2c 13 c 23) — C 13 C 22 


2 (2.97) 


(2-98) 


s 12 — 
s 13 = 


1 

s 44 > 

C 44 

(2.99) 

1 

s 55 — > 

(2.100) 

C55 


1 

s 66 ~ „ - 
c 66 

(2.101) 

“( c 12 c 33 — c 13 c 23) 



c ll( c 22 c 33 — c 23) “ c 12( c12c 33 ~ 2c 13 C 2 3) C ?3 C 22 

c 12 c 23 ~ c 13 c 22 

c ll( c 22 c 33 — c 23) — c 12( c 12 c 33 — 2Ci3C 2 3) — ^\z c 22 


(2.103) 
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Table 2.0 - Calculated linear elastic compliance values using ultrasonic 
stiffness data for transversely isotropic and orthotropic symmetries. 


Compliance Moduli (GPa)" 


Modulus 


Transverse Isotropy 


c, a = 3.0 


c i3 = 8.9 



Orthotropy (c^ = 7.2) 


Modulus 



-0.0049 


-0.0049 


-0.0026 
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and 


s 23 — 


— ( c ll c 23 — c 13 c 12) 


c ll( c 22 c 33 — c 23) — c 12( c 12 c 33 — 2CJ3C23) — Ci 3 C 2 2 


(2.104) 


The calculated values for these moduli are also presented in Table 2.9. Again there was 
a narrow range of values except in the constants s 13 and S23 because of the large uncer- 
tainty in c 13 and C23 . 

The isothermal correction for these adiabatic compliance moduli can be calcu- 
lated using equation (2.34). The temperature of these measurements was nominally 24 
C and the density was 1.54 g/cm 3 . The coefficients of thermal expansion for similar 
unidirectional carbon/epoxy composites were given by Hull [46] as 
q/j = a 2 = 30 ■ 10 -6 /C and a 3 = — 0.2 • 10“ 8 /C and the specific heat was given by as 
c p = 2.1 J/gC (unpublished data from C. Welsh, 1986). Using these values, the max- 
imum correction term was only about 7 • 10“ 6 GPa” 1 which was much smaller than the 
uncertainty of the calculated values. Thus the difference between the adiabatic and the 
isothermal values for the compliance moduli was negligible. 

Some of the elastic compliance moduli were also measured using static compres- 
sive loading with bonded strain gauges to measure the strain. The theory for this is 
based on Hooke’s law. If a unidirectional stress is applied, for example, along the Xj 
axis, the strain can be written as 


6, = Sn(Ti, 

(2.105) 

e 2 = s 12 > 

(2.106) 

e 3 = s 13°l 

(2.107) 


and 


€4=65=66= 0. 


(2.108) 
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Therefore the moduli can be rewritten as 


S 1A 




(2.109) 


or in more general form as 


S AB 


d°A ’ 


( 2 . 110 ) 


for a unidirectional stress along the A axis. Therefore the stress and strain can be meas- 
ured and the slope of the stress-strain curve will provide the moduli. 

The samples used in these measurements were cubes nominally 0.7 in. (1.8 cm.) 
on a side oriented with faces perpendicular to .the coordinate axes. They were also cut 
from the panel MCIS 1 and were designated MCIS 1.8D and E. The load was applied 
with a MTS 810 load frame in compression and the calibrated output of the load cell 
was input into a H.P. 3478A digital multimeter. To insure that the load remained uni- 
directional, a fixture was designed and inserted into the load frame that placed the 
upper compression plate on a ball bearing. This allowed the compression plate to com- 
pensate for any slight lack of parallelism in the specimen. Fig 2.15 illustrates this 
fixture. The strains were measured with Measurements Group Inc. 350Q strain gauges. 
The surfaces of the samples were cleaned according to manufacturer’s instructions and 
the strain gauges were then bonded onto the samples. A conventional Wheatstone 
Bridge circuit was used with a H.P . 6002A power supply providing a five volt input and 
a H.P. 3478A voltmeter measuring the output voltage. A Tektronix 4051 computer was 
used to read the data from the voltmeters and calculate stress, strain and the resulting 
modulus. The control program is shown in Appendix E. Fig 2.16 shows the block 
diagram for the experimental set up. During each test, a unidirectional stress was 
applied and the corresponding strain was measured in two directions. One direction 
was parallel to the direction of applied stress while the other was perpendicular. Then 
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Ball Bearing 

Upper Compression 
Plate 

Specimen 


Lower Compression 
Plate 


Figure 2.15 - Diagram of the compression fixture used during the 
static loading for the measurement of the elastic 
compliance moduli and the uniaxial stress acoustic 
constants. 
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Figure 2.16 - Block diagram of the apparatus used for static 
measurement of the linear elastic compliance 
moduli. 



- 72 - 


the sample was reoriented to apply stress in a different direction. 

The resulting values for the compliance moduli are given in Table 2.10. The shear 
moduli s 44 , S55 and s 6 e were not measured by this technique. Again the deviation from 
transversely isotropic behavior is quite small as demonstrated by the good agreement 
between s u and S22 and between s 13 and S23 for the orthotropic moduli. These values 
compare favorably with the values obtained from the inversion of the ultrasonic data 
except for the value of s 12 . This discrepancy could be caused by error accumulation in 
the inversion of the stiffness moduli of it could be due to inaccuracy of the static meas- 
urement. Ledbetter and Read [47] also compared static measurements to ultrasonic 
measurements for a NbTi fiber reinforced Cu composite and reported discrepancies in 
the off diagonal moduli. The measured values can also be compared to those presented 
by Garber [48] which are listed in Table 2.11. His measurements were made in tension 
on thin laminates of T300/5208. Because thin laminates were used, s 12 was not meas- 
ured. The values compare favorably except for S33 . This is consistent with the lower 
value of C33 reported earlier for the material in the present study. 


II.D. Conclusion 

The linear elastic stiffness and compliance moduli of a unidirectional 
graphite/epoxy composite have been measured by ultrasonic and static measurements. 
These measurements demonstrate a measurable deviation from transversely isotropic 
behavior. However, the deviation is very small, especially when compared to material 
variations from sample to sample, and thus the assumption of transverse isotropy for 
the linear elastic behavior of this material is valid. 

Comparisons of the present data with that presented in the literature for similar 
materials are favorable except for the values of C33 and S33 . The difference in these 
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Table 2.10 - Measured linear elastic compliance values for transversely 
isotropic and orthotropic symmetries. 


Transversely 

Isotropic 

Modulus 

Value 

(GPa) -1 

Orthotropic 

Modulus 

Value 

(GPa)- 1 

s n 

0.089 

s n 

0.0890 

s 12 

-0.063 

s 12 

-0.0626 

s 13 

-0.0021 

s 13 

-0.00208 

s 33 

0.00935 

s 33 

0.00935 



s 22 

0.0891 



s 23 

-0.00215 


Table 2.11 - Comparison of measured compliance values with those 
presented by Garber [49]. 


Modulus 

Present 

Garber 


Work 

(GPa)" 1 

s n 

0.089 

0.0921 

s 12 

-0.063 

not measured 

s 13 

-0.0021 

-0.0025 

s 33 

0.00935 

0.0078 
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moduli together with a lower measured density may indicate a lower fiber volume con- 
tent or a higher void content in the material of the present study. It also may indicate 
either a poor alignment of the fibers or a waviness of the fibers along the fiber axis. 

The comparisons of the ultrasonic and static data are good except for the 
modulus s 12 . The difference may be an artifact of the measurement or the calculations 
used in the comparisons. Because of the higher accuracy of ultrasonic measurements, 
the values used in the computation of the nonlinear moduli in the next chapter are 
those obtained by ultrasonic measurements except in the case of s 13 where large uncer- 
tainties exist. 
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m. NONLINEAR ELASTICTIY 


HLA. Introduction 

Using the linear elastic stiffness and compliance moduli presented in the previous 
chapter, the nonlinear elastic moduli can be calculated from measurements of the stress 
dependence of the ultrasonic "natural” velocity. In this chapter, the theory of how this 
is accomplished will be presented along with actual measurements and calculations for 
the unidirectional graphite/epoxy composites being studied. 


m.B. Theory 

The theory for nonlinear elasticity differs in several respects from the theory of 
linear elasticity which was presented earlier. The first major difference is that the 
assumption of infinitesimal deformation is no longer valid. Therefore, the strain tensor 
must be defined using either the Lagrangian or Eulerian definitions. In this work the 
Lagrangian formulation will be used. The constitutive equation (i.e. relation between 
stress and strain) must be extended to include higher order terms. Also the strain 
energy function must be extended. It is given by Brugger [12] as 

, 1 ,1 

0 2 C ijkl^ij^/kl g c ijklmn^ij^kl^/mn “1“ • (3.1) 


The numerical coefficient for each of the terms is where n is the order of strain for 


n! 


that term. The Voigt reduced notation can be used to rewrite this equation as 


^ ~ o Y c AAVa + £ c ABVaVb + ~Y c AAaVa 


A<B 


+ T £ c AABVaVb + £ c ABD*?A%*7d + • * 


A*B 


A<B<D 


(3.2) 
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where no sum is implied except where indicated and the capital indices range from 1 to 

6. 

The definitions for the elastic moduli are similar to those of linear elastic theory. 
The higher order adiabatic and isothermal moduli were defined by Brugger [12] as 


c ijklmn... Po 


c ijklmn... Po 


s ijklmn... Po 



1 ^ij^'kl^'mn 


(3.3) 

(3.4) 

(3.5) 


and 


T 

& ijklmn... 


= ~ P 0 




[ ^tij^kl^mn J 


(3.0) 


However, in this research, the third order moduli were obtained by measuring the stress 
derivative (isothermal) of the ultrasonic velocity (adiabatic). Therefore the resulting 
moduli were the "mixed" third order moduli as defined by Brugger [12] as 


_ 

''ijklmn 


dcjjki ^ 


dV t 


mn rr 


(3.7) 


or in reduced notation 


C ABD — 




d*k> 


(3.8) 


They are related to the pure adiabatic third order moduli by 


M , T _ S 


(*&] 


c ijklmn * ^ c mnop u: op 

c ijklqr^qr 

<9T 

\ J 

T 


(3.9) 


or 
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g 

C ABD 


nM , 

^ABD + 


PqC'y 


^DE^E 


M 

(*&] 


c ABF Q; F — 

dT 

T 


(3.10) 


However, again the differences are usually small and therefore less than experimental 
uncertainty. 

Using the same arguments based on symmetry of the stress and strain tensors 
and the existence of the strain energy function, the number of independent third order 
elastic moduli can be reduced to 56 for the most general anisotropic material. Further 
considerations of the symmetry of the material reduces the number to 20 for an ortho- 
tropic material and 9 for a transversely isotropic material as shown in Appendix F. For 
the orthotropic case the independent moduli are c lu , c 112 , c n3 , c 122 , 

c 123> C 133» c 144> c 166> C 1M> c 222> c 223> c 233> c 244> c 266> C 2M> c 333> c 344> c 366> c 366 an ^ c 46« * All Other 
moduli are equal to zero. The independent moduli for transverse isotropy are 
c m, c 112 , c n3 , c 123 , c 133 , Cj 44 , c 166 , c^ and c^ with the following relations for the other 
moduli 

c 222 = c lll, c 355 = c 344 

c 223 “ c 1 13, c 166 = c 266 = "^“( c lll “ c 112) 

c 255 = c 144, c 233 = c 133 ( 311 ) 

c 244 = c 155, c 366 = 7r( C 113 "" C 123) 

c 122 “ c 1 12, c 456 5=1 ^"( c 155 ~ c 144) 
with all remaining c = 0 . 

As mentioned previously, these moduli can be obtained from measurements of the 
stress derivative of the ultrasonic "natural” velocity. The theoretical basis for this was 
presented by Thurston and Brugger [14]. They rederived the equation of motion for a 
small amplitude elastic wave propagating in a homogeneously stressed medium using the 
following notation. The coordinates aj refer to the position of a point of material in the 
"natural" or unstressed state. X* are coordinates of the homogeneously stressed material 
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or the material in the initial state while Xj refer to the material in the stressed state 
deformed by the small amplitude elastic wave. The vector components \i { are the com- 
ponents of displacement from the initial state due to the elastic wave and are given by 

Uj = Xj - Xj. . (3.12) 


Uj are the components of displacement referred to the natural state and are related to 
U; by 

'ax,' 


Uj = 


da, 


<U 


iu q . 


(3.13) 


The wave equation was then presented as 


C^Ul 


• • x Is 

A) u j = A J k P m ~^ 


p a m 


(3.14) 


where ” refers to values in the initial state 


and Aj kpm is given by 


s 5XjdX k _ s 

^■jkpm = Oj k t pm + ^ a ^ a . ^pqmi- 


(3.15) 


A solution in terms of the "natural" velocity was chosen of the form 


Uj 


(3.16) 


where N is a unit normal in the natural state. 

This solution implies that a wave with a wave front which is a material plane with unit 
normal N moves from the plane N * “a = 0 to N • ~£ = L 0 in the time L 0 /W . Substitu- 
tion of this solution yielded the propagation conditions 

P 0 W 2 Uj = a/ jk U k (3.17) 


where cjj k is given by 
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^jk “ N r N s (^jkfrs + (^qk ^Vqk)^]Tqs)' (3.18) 

The eigenvector U is normalized so that its components are the direction cosines of the 
particle displacement with respect to the natural state. Differentiating equation (3.17) 
with respect to pressure (p) and evaluating at zero pressure yields 

(a>w 2 ); . o = (u jW ; k u k ) p . o (3.i9) 

which was reduced to the following forms for the conditions of hydrostatic pressure and 
uniaxial compressive stress. For hydrostatic loading, 

— (Po^ 2 )p - 0 = 1 + 2 wFhc + Gjic (3.20) 

where F HC = sJ„,U r U, , (3.21) 

GhC = S MUV c uvprq»NpN q U r U, , (3.22) 

and w = (p 0 W*) 0 = (pv\ = c p s rq ,N p N q U r U, . (3.23) 

For uniaxial loading the equation reduced to 

-(Po W2 )p - 0 = 2 wF uc + G uc , (3.24) 

where F uc and G uc are given by 

F uc = s a T brs M a M b U r U s (3.25) 

and 

^uc = s abuv c uvprqs^aMbNpNqU r U s (3.26) 

and are the components of a unit vector in the direction of stress. 

It is assumed for the uniaxial case that the direction of propagation is always perpen- 
dicular to the direction of stress. 
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These equations provide relationships between the linear elastic stiffness and 
compliance moduli, the third order elastic moduli and the derivative of the "natural" 
velocity with respect to stress. However, the measurements made in this study were of 
the Stress Acoustic Constants (H) which were defined by Heyman [31] and Cantrell [49] 
to be 


AW 

W = (p 0 W 2 )p = o 
A<t (2p 0 W 2 ) 


(3.27) 


Since there were a number of different orientations of stress and propagation and polari- 
zation directions, the different SAC’s were designated with three subscripts as . 
The first subscript gave the direction of applied stress which was limited to one of the 
three axes for these experiments. A zero for the first subscript implies hydrostatic load- 
ing. The second and third subscripts were the directions of propagation and polariza- 
tion of the ultrasonic wave respectively. These again were limited to being pure mode 
waves along the three axes. 

Using this notation, the equations for the SAC’s were derived for both the models 
of transverse isotropy and orthotropy. The results for the 27 cases considered are given 
in table 3.1. 


III.C. Experiment 

The measurements of the stress derivatives of the "natural" velocity were carried 
out using some of the samples that were used in the previously described ultrasonic velo- 
city measurements. Sample MCIS 1.10A was used for the measurements of the longitu- 
dinal velocity derivatives while samples MCIS 1.8A, 1.10B and 1.7A were used in the 


shear wave measurements. 
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Table 3.1 - Stress Acoustic Constant Equations for transversely isotro- 
pic and orthotropic media. 



Transverse Isotropy 


1 4 2 c 11 (s 11 + s 12 + s 13 ) 4 (s n 4 s 12 4 s 13 )(c lu 4 c 112 ) 4 ( 2 s 13 -f 



Ham 


1 + 2055 ( 23 ^ 4* s^) 4* 2(s n 4* S 12 4- s 13 )c 133 + (2s 13 4 




1 4 2c 44 (s n 4 - s 12 4 s 13 ) 4- (s n 4- S 12 4- s 13 )(c 144 4 c 166 ) 4 (2s 13 4 


2^44 




1 -h(cn — c 12 )(s n -j- s 12 4 S 13 ) 4 ^((sn 4 s I2 4 s 13 )(c m — c U2 ) 



1 4 2c 44 (2s 13 4 s^) 4 (s n 4 s 12 4 s 13 )(c 144 4 c 166 ) 4 (2s 13 4 s&)c ^ 4 



2C33?1Z + (sn 4 s 12 )c 133 4 



2c n s 12 4 s u c 112 4 s 12 c in 4 s 13 c 113 




2c u s 13 4 s 13 (c m 4 c 112 ) 4 
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Table 3.1 - Continued 


Transverse Isotropy (Continued) 


H 123 — “ 

2C44 s i3 + s ll c 144 + + s 13 c 344 


2C44 

H 213 — ^123 

h i21 = - 

( c n “ ^12) + l/4(sn 4* S 12 )(c m — c 112 ) 4- 1/2 s 13 (c 113 — C 123 ) 


( C I1 ~ c 12) 

H 212 — H 12 i 

H 312 

( c ll “ c 12) s 13 + ^ s 13( c lll c 1 12 ) + ^ s 33( c 133 “ c 123) 


(c 1 1 — C 12 ) 

H 321 — H 312 

^313 — “ 

2044833 -f S 13 (c 144 -f C 166 ) -f 


2c 44 


H323 — H313 


h 131 = - 

1 

1 

2c 44 s u -f S llCl66 + S 12 C 144 4- S 13 C 344 


2c 44 

H232 — H 1S1 

Him 

2c 44 S 12 4* Sji c i44 -f S 12 C 166 4* S 13 C344 


2 C<4 
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Table 3.1 - Continued 


Orthotropy 

Hon - - 

1 4- 2c 11 S 1 4- SjCm 4- S 2 c 112 4- S 3 Cn 3 


2c n 

Hq22 — — 

1 4- 20^2 4“ SjC^ 4- s 2 c 222 4“ 830^ 


2C22 

Ho 33 = ~ 

1 4- 2C33S3 4* SjC^ 4- + S3C333 


2033 

H ()31 = ~ 

1 4" 2c55Sj 4- SjCj5g 4- S 2 C355 4" S3C355 


2066 

Hq32 ~ 

1 4- 2c 44 S 2 4- SjC 144 4- S 2 c 244 4- S 3 C3 44 


2C55 

Hoi2 = — 

1 4- 2c 66 S 2 4- SjC^e 4- S 2 c 2ee 4- S 3 c 366 


2c«« 


1 4- 4- SjCjgg 4- 4* S 3 C3g fl 


2c« 

^013 = ~ 

1 4" 2C55S3 4- 4- 4- 830335 


2C55 

^023 “ “ 

1 4- 2c 44 S 3 4- S x c 144 4- S 2 c 244 4- S3C355 


2 c 44 

H133 = “ 

- c 33 s 13 4- S U C 123 4- S 12 C233 4* S 13 C 333 
2c 33 


^233 = “ 

2C33 s 23 4* S 12 C 133 4- s 22 c 233 4" s 23 c 333 
2033 


ffi 

s 

II 

1 

2c 2 2 s 12 4- SijC 12 2 4- s 12 c 22 2 4- S 13 C223 
2c 22 


h 2 „ = - 

2CnSi2 4- S 12 C 1U 4- S22C112 4- &23 c \\3 


2c u 


where Sj = s u + s 12 4- s 13 
S 2 = s 12 S22 4 - S03 

and S3 = s 13 4- 4- S33 
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H322 — 


- c 22 s 23 + S13C J22 + s 23 c 222 4 S33C223 




2 c U Si 3 + S 13 C in -1- S23C 112 + S33C 113 


2c u 


2 cq 0 S u 4* s n c 1M 4- s 12 c 2 6 o 4 Si 3 c 360 



H-jio — — 


Hoqi — 


Hi 03 — “ 


H213 — 


c 00 s 23 4 S 13 C 166 + 



c 56 s 23 4 S 12 C 166 + + S25C555 


2c 56 


2C44S55 4* S 13 C H4 4- S23C244 + 8530544 



2c& 6 s 33 + s 13 c 165 + s 23 c 256 4 8550355 



2C6fiSn 4 S U C 156 + S12C255 4* S 13 C355 



2044820 4- Si 2 C 14 4 4- 4- S23C344 



2c 44 s 12 4“ s 11 c 14 4 4- s 12 c 2 44 4- S 13 c 344 




c 56 s 12 4 S 12 C 156 4- 3220 255 + §230555 


2C66 
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The change in "natural" velocity was measured using a pulsed phase locked loop 
(P2L2) interferometer. The theory of the measurement of the change in resonant fre- 
quency and it’s relationship to the change in "natural" velocity was presented earlier. 
Now the theory of the operation of the P2L2 which is shown in Fig. 3.1 will be dis- 
cussed. 

The heart of the P2L2 is a voltage controlled oscillator (VCO) which generates a 
continuous wave signal of a frequency that is controlled by a D.C. input signal. A por- 
tion of the signal from the VCO is gated out into a tone burst and is then used to excite 
the ultrasonic transducer. The resulting elastic wave is launched into the sample and is 
then detected either by the same transducer in the pulse echo mode or another trans- 
ducer in the through transmission mode. The received signal is input into the P2L2 
where it is phase compared with the signal from the VCO at a preselected phase point 
using a sample and hold. The sampled voltage from the phase detector is then used to 
drive the VCO to a condition of quadrature. The acoustic phase shift (0) given by 

6 = 27Tft (3.28) 

is then maintained as a constant by the feedback loop and therefore 

A0 = 0 = 27r(tAf + fAt) (3.29) 

and 

A 6 _ Af At 

e ft' 

Thus, as stated earlier, 

Af = AW = _ At 
f W t ’ 

Therefore, by monitoring the change in the frequency of the P2L2, the change in 
"natural" velocity is also monitored. 


(3.31) 
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phase locked loop 
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A block diagram of the experimental apparatus used for the uniaxial SAC meas- 
urements is shown in Fig. 3.2. The uniaxial compressive load was applied as in the 
strain gauge measurements with a MTS 810 load frame. The calibrated output of the 
load cell was monitored by the H.P. 3478A digital multimeter. The tone burst from the 
P2L2 was used to excite a conventional damped broad band 2.25 MHz transducer. 
Because broad band transducers were used, phase shifts at the transducer specimen 
interface were assumed to be negligible. The P2L2 was used in a pulse echo mode with 
the reflected echos input back into the instrument for phase comparison. A Tektronix 
2445 oscilloscope was used to view the ultrasonic signal, the phase signal and the sample 
and hold position. It was triggered by a reference sync from the P2L2. The frequency 
was monitored by a H.P . 5316A universal counter. The voltage and frequency values 
were read by a Tektronix 4051 computer and converted to stress and normalized fre- 
quency shift. The program listed in Appendix G also computed the slope of the curve 
by linear regression as well as plotted and stored the data. 

The value for each uniaxial SAC was determined by averaging the values from 
ten measurements. A spring loaded clamping device was used to maintain constant 
pressure of the transducer on the specimen in order to hold the bond thickness constant 
during the measurement. Variations in the bond thickness would cause phase shifts 
which introduced errors into the measurement. These errors which have less effect on 
thick samples made initial measurements on thin laminates impossible and placed 
another constraint on minimum specimen thickness. The transducer was also rebonded 
onto the specimen three or four times during the ten measurements to determine the 
effect of slight variations in initial bond thicknesses. After each rebonding, the specimen 
was ramped through the load cycle several times to settle the bond. Then, the sample 
was allowed to remain at rest approximately 15 minutes to bring the sample to equili- 
brium temperature. This was necessary because the handling of the sample during 
rebonding increased the temperature slightly which could significantly affect the results. 


C -3- 
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Figure 3.2 - Block diagram of uniaxial stress acoustic constant 
measurement apparatus. 
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Another concern in these measurements was again the presence of viscoelastic 
effects in these materials. The loading rate was varied considerably to examine this 
effect. The specimen were loaded to between 80 and 100 MPa in times ranging from 10 
sec. to over three minutes with no observed effects. The measurements were made for 
both increasing and decreasing load with no measurable differences. The effect of tem- 
perature variations due to thermoelastic heating was also considered. Theory assumes 
these measurements to be made under isothermal conditions. At these loading rates 
which were very nearly isothermal, the effect of thermoelastic heating was calculated to 
be much less than experimental uncertainty. The basis of these calculations was the 
temperature derivative of velocity data to be presented later. 

The previously mentioned constraints on the choices of stress, propagation and 
polarization directions allows for 18 measurements to be made. The average values 
from the ten measurements of each and their standard deviations are presented in 
Table 3.2. The uncertainties were typically less than five percent of the measured 
values for these measurements. The values were consistently positive indicating increas- 
ing frequency as a function of increasing compressive stress except for H 323 and H 3 i 3 . A 
few of the wave modes exhibited nonlinear (quadratic) curves indicating higher than 
third order nonlinear effects. These are designated by an * in the table and the SAC 
values for these are the linear coefficients of a quadratic fit to the data. Typical curves 
for the uniaxial SAC measurements are shown in Figs. 3.3 - 3.20. 

As with the linear elastic measurements, some of these measurements can be used 
to check for consistency with transversely isotropic behavior. The SAC pairs 
(H122, Ha,,) (H 311 , II322) , (H 133 , H233) , (H„„ H232) , and (H 132 , H^i) should all have equal 
values. Other pairs that should have equal values were not compared as they were 
measured on different specimens making comparisons invalid because of sample to sam- 
ple variations. The values of some of these pairs compare well within experimental 
uncertainty. However, H 133 and differ by more than the uncertainty of the 
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Uniaxial Stress Acoustic Constants 


Value (GPa)" 


0.0490 +/- 0.0009 


0.0427 +/- 0.0010 


0.00116 +/- 0.00007 


0.00123 +/- 0.00001 


0.0887 +/- 0.0007 


0.0741 +/- 0 001 ° 


0.00279 +/- 0.00010 


0.00299 +/- 0.00010 


0.068 +/- 0.002 


0.0572 +/- 0.0010 


- 0.00993 +/- 0.0004 


- 0.00919 +/- 0.0003 


0.165 +/- 0.003 


0.149 +/- 0.002 


0.111 +/* 0.002 


0.109 +/- 0.001 


0.0538 +/- 0.0005 


0.0479 +/- 0.0003 
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Figure 3.3 - Plot of normalized frequency change as a function of 
uniaxial stress for a longitudinal wave propagating 
along Xj with stress applied along X3. 



STRESS (Pa) 

Figure 3.4 - Delta F/F versus stress for a longitudinal wave 

propagating along x 2 with stress applied along x-j. 
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Delta F/F versus stress for a longitudinal wave 
propagating along X 3 with stress applied along x 
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Figure 3.6 - Delta F/F versus stress for a longitudinal wave 
propagating along X 3 with stress applied along x 
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Figure 3.7 - 


Delta F/F versus stress for a longitudinal wave 
propagating along X£ with stress applied along 
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Figure 3.8 - 


Delta F/F versus stress for a longitudinal wave 
propagating along with stress applied along x 
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Figure 3.9 - Delta F/F versus stress for a shear wave propagating 
along x 2 polarized along x^ with stress applied along 


X l* 



Figure 3.10 - Delta F/F versus stress for a shear wave propagating 
along x^ polarized along X£ with stress applied along 
x 2 . 
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3,11 - Delta F/F versus stress for a shear wave propagating 
along x^ polarized along X 2 with stress applied along 

x 3* 



Figure 3,12 - Delta F/F versus stress for a shear wave propagating 
along X 2 polarized along x^ with stress applied along 
x^. 
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Figure 3.13 - Delta F/F versus stress for a shear wave propagating 
along x 2 polarized along x^ with stress applied along 


X l* 



x 2 . 
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Figure 3.15 - Delta F/F versus stress for a shear wave propagating 
along X 2 polarized along x^ with stress applied along 
x^. 
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Figure 3.16 - Delta F/F versus stress for a shear wave propagating 
along Xj polarized along x^ with stress applied along 
x r 
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Figure 3.17 - Delta F/F versus stress for a shear wave propagating 
along polarized along x^ with stress applied along 

x l* 



Figure 3.18 - Delta F/F versus stress for a shear wave propagating 
along x 3 polarized along x 2 with stress applied along 
X£ . 
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Figure 3.19 - Delta F/F versus stress for a shear wave propagating 
along polarized along x 2 with stress applied along 

X l* 



Figure 3.20 — Delta F/F versus stress for a shear wave propagating 
along X 3 polarized along Xj with stress applied along 


x 2 . 


- 100 - 


measurement as well as H 131 and H232 . Thus, again measurable deviations from 
transverse isotropy were demonstrated. However, why this behavior is demonstrated 
only in some measurements and not in others is not clear. 

The hydrostatic measurements were made in a similar fashion as those made 
under uniaxial compression. A block diagram of the set up is shown in Fig. 3.21. The 
change in "natural" velocity was again measured with the P2L2. Argon was used to 
increase the pressure within the pressure chamber. A HP. 5316A universal counter 
monitored the frequency while a Tektronix 465B oscilloscope was used to view the ultra- 
sonic, phase, and sample/hold signals. The Tektronix 4051 computer was used to read 
the frequency data. The pressure data was manually read from a Heise C-57488 pres- 
sure gauge and input into the computer. 

There were several differences from the uniaxial measurements, however. While 
temperature variations during the uniaxial SAC measurements were insignificant, they 
were important in the hydrostatic measurements. Increasing the pressure 25 psi. sud- 
denly, caused an increase in temperature of the gas and thus the sample of several 
tenths of a degree C. This temperature increase was significant. Therefore, the pres- 
sure had to be increased slowly and then time for thermal equilibrium to be obtained 
had to be given. Typically, there was a 20 to 30 minute delay between data points. 
Also, over this period of time the room temperature could vary significantly. This 
caused additional problems. To overcome this, heating elements were installed in the 
pressure chamber. A temperature controller was used with these heaters to maintain 
the temperature several degrees C above room temperature. This speeded the time 
required to attain thermal equilibrium and avoided the effect of room temperature vari- 
ations provided it did not rise more than the several degree buffer. A Lake Shore Cryo- 
tronics Inc. model DRC84C temperature controller was used. This instrument had a 
reported ability to maintain the temperature within +/- 0.1 degree C. Although this 
instrument drastically improved the situation, these small variations in temperature 
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Figure 3.21 - Block diagram of apparatus used in the hydrostatic 
pressure stress acoustic constant measurements. 
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reduced the accuracy of the measurement significantly. 

Another problem that increased the uncertainty in the hydrostatic measurements 
was the limitation in pressure. The upper limit of pressure in the chamber was set at 
only 250 psi. This was about a factor of 50 less than the maximum stress applied dur- 
ing the uniaxial measurements. Theoretically, the slope of the SAC should be evaluated 
at zero stress and thus this limitation in pressure should not be a problem. However, the 
low pressure produced much smaller velocity changes making them more uncertain. 
Also, the effect of temperature variations were much greater. Thus, the uncertainty of 
these measurements was much larger. It was on the order of ten percent of the meas- 
ured SAC for most of the hydrostatic measurements. The limitation in maximum pres- 
sure also made the measurement of impossible because the changes involved were to 
small to detect with the P2L2. The possibility of anomalous behavior at low pressures 
must also be considered. Therefore, although the hydrostatically measured SAC values 
were used with the uniaxial SAC measurements to determine the nonlinear moduli in 
this work, measurements at higher pressures are needed to determine the validity of 
these calculations. 

Because of the tremendous amount of time involved in making a single hydros- 
tatic SAC measurement (10 to 15 data points at 20 to 30 minutes per data point), only 
three measurements were used in the calculation of each hydrostatic SAC. This also 
contributed to the large uncertainty in these measurements. Measurements were again 
made for increasing and decreasing pressures with no observed differences. 

The hydrostatic measurements were all made using undamped crystal transduc- 
ers (PZT 5A) because of size limitations within the pressure chamber. This created the 
necessity of making corrections because of phase shifts at the transducer-specimen inter- 
face. These phase shifts are a function of frequency and can be calculated using the 
theory presented by McSkimmin [19]. The correction to the experimentally determined 
normalized frequency shift assuming zero bond thickness was given by 
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where F 0 , F l were the measured frequencies of two data points, 
to was the time of flight in the unstressed specimen, 
and 7 0 , ^ were the phase angles of the two data points in degrees. 
They were given by 
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(3.35) 


where Z ( , Z x were the mechanical impedances of the specimen and 
transducer respectively, 

and 


Rq — Zy tan 


180 


'R 


(3.36) 


and 
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R t = Z x tan 



(3.37) 


where Fr was the resonant frequency of the transducer. 

The effect of nonzero bond thickness on this correction term was also analyzed and 
found to be less than experimental uncertainty. 

The original and corrected values for the hydrostatic SAC’s are shown in Table 
3.3. Typical curves are shown in Figs. 3.22 - 3.29. Comparisons of the data to check for 
transversely isotropic behavior show no measurable differences. This is most likely due 
to the large uncertainty of the measurements. 

The effect of temperature on the "natural" velocity was also measured for small 
variations about room temperature. The apparatus used for these measurements was 
the same as that used in the hydrostatic measurements. The temperature was increased 
above room temperature several degrees C in 0.3 C steps using the heaters and the tem- 
perature controller. The frequency was measured at each point after thermal equili- 
brium was attained. Again corrections were applied to account for phase shifts at the 
transducer-specimen interface. Typical curves are shown in Figs. 3.30 - 3.38. The slopes 
of these curves were designated as the Thermal Acoustic Constants (Hxy) and are given 
in Table 3.4 where the second and third subscripts give the directions of propagation 
and polarization respectively. The same samples used for the hydrostatic and uniaxial 
measurements were also used for these measurements. 

The curves all demonstrated linear relationships between temperature change 
and normalized frequency shift for the limited range of the measurements. While these 
measurements were not used in the data reduction to determine the third order elastic 
coefficients, they were important in determining the effect of temperature changes on 


the SAC measurements. 
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Table 3.3 - Measured and corrected values for the Hydrostatic Stress 
Acoustic Constants. 


Hydrostatic Stress Acoustic Constants (GPa) -1 

SAC 

Original Value 

Corrected Value 

Hon 

0.46 +/- 0.02 

0.39 +/- 0.02 

Hq22 

0.44 +/- 0.02 

0.37 +/- 0.02 

H(H2 

0.301 +/- 0.007 

0.257 +/- 0.007 

Hq21 

0.305 +/- 0.006 

0.261 +/- 0.006 

Hoi3 

0.32 +/- 0.03 

0.27 +/- 0.03 

H023 

0.32 +/- 0.02 

0.27 +/- 0.02 

H<m 

0.33 +/- 0.02 

0.28 +/- 0.02 

H()32 

0.30 +/- 0.02 

0.25 +/- 0.02 

Hq33 

Not Measured 
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Figure 3.22 - Delta F/F versus pressure for a longitudinal wave 
propagating along x^. 



STRESS (Pa) 

Figure 3.23 - Delta F/F versus pressure for a longitudinal wave 
propagating along x 2 » 
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Figure 3.24 - Delta F/F versus pressure for a shear wave propa- 
gating along polarized along x 2 . 



Figure 3.25 - Delta F/F versus pressure for a shear wave propa- 
gating along x 2 polarized along x^. 
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Figure 3.26 - Delta F/F versus pressure for a shear wave propa- 
gating along xj polarized along x^. 



Figure 3.27 - Delta F/F versus pressure for a shear wave propa 
gating along x 2 polarized along X3. 


AF/ 
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Figure 3.28 - Delta F/F versus pressure for a shear wave propa- 
gating along polarized along Xj, 



Figure 3.29 - Delta F/F versus pressure for a shear wave propa 
gating along x^ polarized along X2« 
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Table 3.4 - Measured and corrected values for the Thermal Acoustic 
Constants. 


TAC 


H X n 


H-T22 


H T 12 


Ht21 


Hti3 


Hj23 


Ht31 


h X32 


Ht33 


Thermal Acoustic Constants (10 -4 


Original Value 


-7.58 +/- 0.02 


-7.55 +/- 0.03 


-9.9 +/- 0.1 


-9.9 +/- 0.1 


-11.9 +/* 0.2 


- 11.8 +/- 0.2 


-8.4 +/- 0.5 


-9.0 +/- 0.7 


-0.84 +/- 0.09 


c- 1 ) 


Corrected Value 


-6.36 +/- 0.02 


-6.33 +/- 0.03 


-8.5 +/- 0.1 


-8.5 +/- 0.1 


- 10.2 +/- 0.2 


- 10.1 +/- 0.2 


-7.1 +/- 0.5 


-7.6 +/- 0.7 


-0.72 +/- 0.09 
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3.30 - Delta F/F versus temperature change for a longitudinal 
wave propagating along Xp 
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Figure 3.31 - Delta F/F versus temperature change for a longitudinal 
wave propagating along X2« 
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3.32 - Delta F/F versus temperature change for a shear 
wave propagating along x^ polarized along x^. 



Figure 3.33 - Delta F/F versus temperature change for a shear 
wave propagating along x^ polarized along X 2 * 
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3„34 - Delta F/F versus temperature change for a shear 
wave propagating along x^ polarized along X£. 
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Figure 3.35 - Delta F/F versus temperature change for a shear 
wave propagating along X 2 polarized along x^. 
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re 3.36 - Delta F/F versus temperature change for a shear 
wave propagating along x^ polarized along X 3 . 
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Figure 3.37 - Delta F/F versus temperature change for a shear 
wave propagating along X 2 polarized along x^. 
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Figure 3.38 - Delta F/F versus temperature change for a longitudinal 
wave propagating along x^. 
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m.D. Data Analysis 

The SAC data, along with the linear moduli were used to determine some of the 
third order elastic moduli for this material. Since there were more equations and data 
than unknown moduli, the data was reduced by using a least squares fit to determine 
the best values for the moduli. The data was fit using both the transversely isotropic 
and orthotropic models. However, the lack of data for prevented the calculation of 
the moduli c 133 and C 333 for the transversely isotropic model and the moduli 
c i33> c 233> c 333 > c i23> an d C 4W f° r orthotropic model. 

The least squares procedure for the determination of third order moduli was out- 
lined by Hankey and Schuele [40], The formulation is as follows. The measured values 
for the SAC’s are added to the constant terms in the equations of Table 3.1. These are 
designated Mj . The coefficients of the equations are designated Ay and the moduli to 
be determined are Dj. Thus the equations of Table 3.1 are rewritten as 

Mi = AjjDj (3.38) 

where i is summed over the number of data used and j is summed over the number of 
moduli to be determined. This can be rewritten as 


A lk Mi = (A ik A ij )D j (3.39) 

where k also runs from one to the number of coefficients to be determined. The best fit 
for the moduli Dj are then given by 

Dj = (Aj k Ajj) A pk Mp (3.40) 

where (A^Ay) -1 is the inverse matrix element and p is summed over the number of 
measured values. 

If it is assumed that the error in the measured values is much greater than the 
error in the coefficients, then the error in the coefficients can be dropped. The error in 
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the fitted values is then given by 


ADj - KAjA.r'ApklAMp 


where AD, are the maximum limits of error, 
and AM p are the errors in the measured values. 
The probable limit of error (ADj) is given by 



ADj 

(N - 1)* 


(3.41) 


(3.42) 


where N is the number of measurements. 

The above stated assumption was valid in these measurements as the errors in the 
measured SAC’s ranged from one to ten percent while the errors in the linear coefficients 
were typically less (parts in 10 3 ) 

Programs were written to perform the data reduction for the two symmetry 
models and are listed in Appendix H. The fitted values and their uncertainties are 
given in Table 3.5. 


HI.E. Conclusion 

Using a pulsed phase locked loop interferometer, the normalized change in the 
ultrasonic "natural” phase velocity as a function of both stress and temperature have 
been measured. The stress measurements were made under conditions of hydrostatic 
and uniaxial stress. The measured curves were predominantly linear although several 
were nonlinear (quadratic) indicating higher than third order nonlinear elastic effects. 
Because of equipment limitations, the hydrostatic measurements were made at much 
lower pressures in comparison to the uniaxial measurements. Thus the uncertainties of 
the hydrostatic measurements were much larger making the validity of using them with 
the uniaxial measurements to calculate the third order moduli questionable. This 
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Table 3.5 - Calculated values of third order elastic stiffness moduli for 
transversely isotropic and orthotropic symmetries. 
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limitation in maximum hydrostatic pressure also prevented the measurement of 
which prevented the calculation of several of the nonlinear moduli. The hydrostatic 
measurements need to be evaluated at higher pressures to reduce the uncertainties and 
to evaluate H033 . This would enable the determination of the remaining third order 
moduli. 

The effect of temperature changes on "natural" velocity was also evaluated for 
nine combinations of propagation and polarization directions. These measurements 
were made over a narrow range of temperature change slightly above room tempera- 
ture. The curves all demonstrated linear relationships between velocity and tempera- 
ture. Although not directly used for the evaluation of third order moduli, these meas- 
urements were important in the evaluation of the effects of temperature variations dur- 
ing the SAC measurements. 

The SAC values and the second order moduli were used to calculate some of the 
third order elastic moduli. A least squares fit was applied for both the transversely iso- 
tropic and the orthotropic models of elastic symmetry. The values for these moduli 
were all negative with the exception of the C123 • This is in agreement with the non- 
linear behavior of conventional materials. Although the linear elastic moduli and the 
SAC data demonstrated small deviations from transverse isotropy, large uncertainties in 
the third order moduli masked any such deviations. The orthotropic moduli agree 
within experimental error to the conditions of transverse isotropy. 
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IV. SUMMARY AND CONCLUSIONS 


The mechanical behavior of a unidirectional T 300/5208 graphite/epoxy compo- 
site material has been evaluated in this research. In particular, the linear (second 
order) and nonlinear (third order) elastic properties have been measured. The linear 
moduli were measured by ultrasonic velocity measurements as well as by static loading 
strain gauge measurements. The nonlinear moduli were calculated from the linear 
moduli and measurements of the change in ultrasonic "natural" phase velocity as a 
function of stress. 

The measured linear elastic properties of this material agreed well with previ- 
ously reported values for a similar material except for the moduli C33 and S33 . These 
were the moduli that reveal the behavior of the composite under loading along the fiber 
direction. The low value for C33 and corresponding high value for S33 may have been a 
result of a lower fiber content or higher void content in this material. It could also have 
been the effect of misalignment of fibers in some of the laminae. However, the values of 
the ultrasonic and static measurements compared well with the exception of the off 
diagonal modulus c 13 . This discrepancy was probably the result of either large uncer- 
tainties of the static measurements or error propagation in the inversion of the ultra- 
sonic moduli. These measurements also pointed out measurable deviations from 
transversely isotropic elastic behavior which is usually assumed for unidirectional com- 
posites. The material more correctly possessed orthotropic symmetry although the devi- 
ations from transverse isotropy were small. 

The change in ultrasonic natural' velocity as a function of stress and tempera- 
ture was measured with a pulsed phase locked loop ultrasonic interferometer for 
numerous wave modes. These measurements again demonstrated deviations from 
transversely isotropic behavior. Also, some of the stress- velocity curves were nonlinear 
(quadratic) which indicated higher than third order nonlinear effects. The 
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temperature-velocity curves were all linear over the small range of the measurements. 
Some of the third order elastic stiffness moduli were calculated from the stress depen- 
dence of velocity. These were all negative in agreement with nonlinear properties of con- 
ventional materials. The uncertainties of some of these moduli were large as a result of 
the large uncertainties of the hydrostatic pressure measurements propagating in the 
least squares data reduction. These large uncertainties masked any deviations from 
transverse isotropy that might have been present. 

However, it should be noted here that the measured SAC’s are probably more use- 
ful than the calculated nonlinear moduli. This is because they are more easily obtained 
with smaller experimental uncertainty and can be measured in an individual specimen 
instead of calculated from measurements made in several specimens. Also, the SAC’s are 
the quantities that are more useful in experimental attempts to measure applied and 
residual stress and may be important in nondestructive evaluation of ultimate strength 
of materials. 

It would be remiss to fail to point out several factors that could be improved to 
yield better measurements. First, the hydrostatic measurements should be carried out 
to higher pressures with improved temperature controller capabilities. This would 
reduce the uncertainties of these measurements and thus those of the nonlinear moduli. 
Also, the effect of sample to sample variations on both the linear and nonlinear proper- 
ties needs to be evaluated further. These measurements were made on as few samples 
as possible to reduce this effect in the calculations. Future measurements should also 
evaluate the changes in velocity over wider temperature variations and with much 
higher stresses (i.e. near ultimate strength stresses) as well as evaluate the nonlinear 
properties by other techniques such as harmonic generation. 

Although the primary goal of this investigation was a more thorough characteri- 
zation of the mechanical properties of graphite/epoxy composites, these measurements 
serve as a basis for further study and potential applications in the area of 
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nondestructive testing of these materials. Also, the effects of fiber volume ratio and 
more complex lay up orientations on nonlinear properties of composites need to be 
evaluated. The relationships between nonlinear properties and important engineering 
properties such as strength, residual strength after impact and fatigue loading, and 
fiber-matrix interfacial strength need to be studied. These may lead to the development 
of needed nondestructive evaluation techniques for composite materials which will play 
an important role in the future of the aerospace industry. 

The problem or residual strain in composites which arises because of the 
mismatch of coefficients of thermal expansion between fiber and matrix is another 
potential application for these measurements. When the composite material is cooled 
from the higher temperatures of polymerization during cure, large strains develop which 
may degrade the ultimate strength of the material. Although special curing cycles have 
been developed to reduce this effect, it still remains as an important problem in the field 
of composites. And as nonlinear ultrasonics have been used in attempts to nondestruc- 
tive^ evaluate residual strains in conventional materials, they may well be useful in the 
measurement of residual stains in composites. Especially, since it may be possible to 
characterize the elastic symmetry of the composite better than in conventional materi- 
als because the symmetry of the material is determined by the lay up orientation. In 
conventional materials the problems associated with anisotropic textures of the material 
have limited the applicability of this technique. 
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appendex a 


The derivation of the independent linear (second order) elastic moduli for the 
symmetries of transverse isotropy and orthotropy are now presented. The case of ortho- 
tropic elastic symmetry is presented first. In an orthotropic material, there are three 
orthogonal planes of two fold symmetry. Therefore if the three coordinate axes are 
chosen such that they are the normals to these planes, then the independent elastic 
moduli can be determined by applying the two fold rotations about each of these axes 
on the strain energy function. The strain energy function is initially given by 

$ = — (cjjcf + c 22 e| + C 33 C 3 + c 44 e| + c 55 e 5 2 + c 66 e|) + 

2 

+ C 13 c 1 6 3 + C 14 € 1 € 4 + + C 23 e 2 e 3 + C 24 e 2 e 4 + c 25 e 2 e 5 

+ C 26 e 2 e 6 + c 34 e 3 e 4 + c 35 e 3 e 5 + c 36 e 3 e 6 + c 45 e 4 e 5 + c 46 e 4 e 6 + c 56 e 5 e 6 . (A.i) 


The effect of a 180 degree rotation about the x 2 axis is now considered. The 
matrix of transformation is given by 



-10 0 
0 10 . 

0 0-1 


(A.2) 


The effect on this rotation on the strains is given by 

/ 

e ij “ a im a jn e mn 


(A.3) 


which yields the following relations between the rotated and original strains expressed 
in the Voigt reduced notation 

4 = C 2 
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The strain energy function can be written for the rotated state and the above equations 
substituted to yield 

= — (c^ 2 + C 22 e| + c 33 e 3 + c 44 e 4 + c 55 e 5 + c 66 e l) + c 12 e l e 2 

2 

+ c 13 € 1 € 3 — c 14 e l e 4 + c 15 e l C 5 — c 16 e l e 6 + c 23 € 2 e 3 “ c 24 e 2 e 4 + c 25 e 2 e 5 

“ C 26 «2 e 6 ~ c 34 € 3 € 4 + c 35 c 3 e 5 “ c 36 e 3 e 6 “ c 45 e 4 e S + c 46 e 4 e 6 “ c S6 e 5 e 6- ( A 5 ) 

Equating <f> and 4>' term by term with respect to strain yields the following 

c 14 “ — c 14 
c 16 = c 16 
c 24 = — c 24 

C 2 6 = — c 26 

. (A.6) 

c 34 — — c 34 
c 36 ~ — c 36 
c 45 = C 45 

c 56 = c 56 

For these relations to be true, then 

c 14 = C 16 = C 24 = c 26 = c 34 = c 36 = c 45 = c 56 = 0. (A.7 

The matrix of elastic stiffness moduli is then given by 
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[ c ab1 


C 11 

c 12 

c 13 

0 

c 15 

0 

c 12 

c 22 

c 23 

0 

c 25 

0 

c 13 

c 23 

c 33 

0 

c 35 

0 

0 

0 

0 

c 44 

0 

c 46 

c 15 

c 25 

c 35 

0 

c 55 

0 

0 

0 

0 

c 46 

0 

c 66 


(A.7) 


The strain energy function is now given by 

^ = 7^( c ll c l 2 + c 22 e 2 + c 33 + e 3 + e 4 + c 55 e 5 + e !) + c 12 e l e 2 


+ c 13 e l e 3 + c 15 e l e 5 + c 23 e 2 e 3 + c 25 e 2 e 5 + c 35 e 3 e 5 + c 46 e 4 e 6- (A.8) 


The effect of a 180 degree rotation about the axis can now be considered. The 
transformation matrix is given by 



1 

0 

0 


0 0 

-1 0 

0 -1 


(A.9) 


The rotated and original strains are related by 


/ 



(A.10) 


Again these are substituted into the equation for the strain energy function resulting in 


^ ~ 2 ( Clie ! 2 + C 22 e 2 2 + c 33 e ! + c 44 e 4 + c 55 e 5 2 + c 66 e f) + c 12 e l e 2 
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+ c 13 C l e 3 — c 15 e l e 5 + C 23 C 2 C 3 — c 25 e 2 e 5 — c 35 e 3 c 5 — c 46 e 4 € 6- 


Equating term by term yields 

c 15 = c 25 “ c 35 = c 46 = 0 
and reduces the elastic stiffness matrix to 


C 11 c 12 c 13 0 0 0 

c 12 c 22 c 23 0 0 0 

c 13 c 23 c 33 0 0 0 

0 0 0 c 44 0 0 

0 0 0 0 c 55 0 

0 0 0 0 0 c 66 


(A.13) 


(All) 


(A.12) 


The third symmetry rotation about x 3 can be shown to have no further effect on the 
elastic stiffness tensor and thus the final form for an orthotropic body is given by equa- 
tion A.13. There are nine independent moduli for an orthtropic material. 

The reduction of moduli for a transversely isotropic material can be carried out 
beginning with the matrix form for orthotropy. This is because it possesses the same 
symmetry conditions as an orthotropic material in addition to having an isotropic plane 
perpendicular to x 3 . Thus the strain energy function can be written out with the ortho- 
tropic moduli and the effect of an arbitrary rotation around the x 3 axis can be con- 
sidered. $ is given by 

^ = + c 22 e 2 + c 33 e 3 + c 44 e 4 + c 5S e 5 + c 66 e !) 

+ c i2 e i c 2 + c 13 e 1 e 3 + c 23 e 2 e 3 . ( A14 ) 


The effect of an arbitrary rotation around x 3 by an angle 0 on the strains is given by 
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e[ = cos 2 0 e l + sin0cos0e 6 + sin 2 0e 2 
e 2 = sin 2 ^! + sin0cos#e 6 + cos 2 0 e 2 


e 3 = e 3 


c 4 = — sin^Cg + cos$c 4 


(A. 15) 


€5 = cos#e 5 — sin0e 4 

= — 2sin0cos0e 1 + (cos 2 0 — sin 2 0)e 6 + 2sin0cos0e 2 . 


These equations can be substituted into the strain energy function as before. Then com- 


paring $ and term by term will yield 

C 11 — c 22 
c 13 — c 23 
c 44 = c 55 

c 66 " -^(Cll “ c 12) 


(A.16) 


which reduces the number of independent elastic moduli to five for a transversely isotro- 
pic body. The elastic stiffness matrix can be written out as 


C 11 c 12 c 13 0 0 0 

c 12 C 11 c 13 0 0 0 

c 13 c 13 c 33 0 0 0 

0 C 44 0 0 

0 0 C 44 0 


( c ab) = 


0 0 

0 0 


0000 


0 



(A.17) 
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APPENDIX B 


Using the equations derived for linear elastic wave propagation presented in 
Chapter IL, the expressions for the velocities of ultrasonic waves in terms of elastic 
moduli and density are now derived. The expressions are derived for both transverse 
isotropy and orthotropy. The case of transversely isotropic elastic symmetry is con- 
sidered first. 

For the case of transverse isotropy, equation (2.84) can be used to show that 
there are several directions of propagation for pure mode elastic waves. Pure mode 
elastic waves may be propagated along the x 3 axis and along any direction perpendicu- 
lar to the x 3 axis in the isotropic plane which of course includes the x t and x 2 axes. 
They may also be propagated in a direction which makes an angle © with the x 3 axis 
given by 


tan 0 = 



+ 2c 44 - c 33 




Cl 3 + 2c 44 - c n J 


(B.l) 


However, for unidirectional graphite epoxy composites where the modulus c^ is usually 
about two orders of magnitude larger than the other moduli, this direction does not 
exist. Thus the pure mode directions to be considered are along each of the three coor- 
dinate axes. 

The components of the matrix X ik for a transversely isotropic material are given 
by 

\l = c lllf + — c 12)l| + c 44^3> 

^22 = ~ Cl2 ^! 2 + + c 44 * 3 > 


^33 — c 44^1 + c 44^2 c 33 l 3 > 


(B.2) 
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Xi2 — ^21 ~ g C 12)Ul2> 


^13 ~ ^31 ~ ( c 13 + c 44)h^3> 


^23 — ^32 ~ ( c 44 + c 13)12*3- 


The solutions for the wave modes propagating along the x t axis are now derived. For 
waves propagating in this direction, the direction cosines of the wave normal are given 
by lj = 1, and 1 2 = 1 3 = 0 . The components of the Xy matrix reduce to 


\i — c n> 


X22 — 2 ( c u “ c 12)* 


^33 — c 44> 


and all other Xy = 0 . These are substituted into equation (2.55) to determine the solu- 
tions for the three modes of propagation in this direction by solving the equation 

C 11 - Po y2 0 0 

DET 0 j(c„ - c 12 ) - p 0 v 2 0 =0. (B.4) 

0 0 C 44 p 0 v 2 


This can be expanded to 


(Cll - Po v2 )(7r( c n - C 12) - Po v2 )( c 44 - Po y2 } = 0. 


The solutions are then found to be 


Cll 2 
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v 2 


C 11 — c 12 




( 2p 0 


(B.6) 


and 


v 3 “ 



ri 


Po ) 


Solving for the eigenvectors for each of these eigenvalues yields the following particle 
displacements. For the first wave, the direction cosines of the particle displacements are 
a 2 = a 3 — ^ and = 1 . Therefore, since = li , the first wave is a pure mode longitu- 
dinal wave. For the second wave, the direction cosines are a x = a 3 = 0 and a 2 = 1 . 
Thus, since = 0 , it is a pure mode transverse wave. The particle displacement vec- 
tor direction cosines for the third wave are found to be ct x = a 2 — 0 and a 3 = 1 , which 
means it too is a pure mode shear wave. The energy flux equations can be used to show 
that there is no energy flux deviation for these waves. 

The solutions for propagation along the x 2 axis are similar. The direction cosines 
of the wave normal are lj = I3 = 0, and 1 2 = 1 . The Xy components are then 

\l = "^"( C ll — C 12)> 

^22 = c ll> (B.6) 

^33 “ c 44> 

and all other Xy — 0 . The solutions for these wave modes are of course the same as for 
propagation along X! except the first and second waves are reversed. Waves propagat- 
ing in any direction in the x l7 x 2 plane with I3 = 0 can be shown to have similar solu- 
tions to the waves propagating along x { and x 2 . 

Propagation along x 3 where li = 1 2 = 0, and 1 3 = 1 is now considered. The com- 


ponents of Xy are now given by 
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\l = c 44> 
^22 = c 44> 
^33 = c 33> 


and all other Xy = 0 . Substitution into equation (2.55) yields 

(C 4 4 - P0 y2 ) 2 ( c 33 “ P0 y2 ) = 0. 


(B.7) 


(B.8) 


The solutions for this equation are given by 


Vl = v 2 = 


c 44 


K Po , 


(B.9) 


and, 


Vo = 


c 33 

Po 


Solutions for the particle displacements for these wave modes show that for the first and 
second wave modes, a 3 = 0 while and a 2 are arbitrary as long as af -f of = 1 • Thus 
these waves are pure mode shear waves with particle displacements in the Xj, x 2 plane. 
For the third wave, they are Q l = a 2 = 0, and 0f 3 = 1 , which means it is a pure mode 
compressional wave. Again, it can be shown that these waves suffer no energy flux devi- 
ation. 

Measurements of the velocities of the waves in these three directions can be used 
to determine four of the five independent elastic stiffness moduli. However, the modulus 
c 13 can not be determined. To evaluate this modulus, an off axis non-pure mode wave 
must be used. The wave chosen for this purpose was a wave propagating in the x lf x 3 
plane 45 degrees between the two axes. Thus, ^ = 1 3 = and 1 2 = 0 . The com- 

ponents of the Xjj matrix are then 
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X 11 = ^"( C ll + C 44)» 

X 22 = ^( c 44 + y( c ll “ c 12))> 

X 33 * “(C33 + C44), (B.10) 

X 13 = "^"( c 13 + C44)' 


and all other Ay = 0 . These are substituted into the eigenvalue equation and the solu- 
tions are given by 


Po v l — 4 [ c ll + 2c 44 + c 33 + \/( c ll — c 33) 2 + 4 ( c 13 + c 44) 2 . 


P 0 V 2 


2 _ _L 


c 44 + “(<>11 — C 12) 


(B.ll) 


and 


Po v I — 4 [ c n + 2c 44 + c 33 — V( c n ~ c 33) 2 + 4(c 13 + c 44 ) 2 . 


Solutions for the particle displacements of these waves show that the second wave is a 
pure mode shear wave. However, the first wave is a quasilongitudinal wave and the 
third wave is quasishear. Energy flux equations show that these wave modes also suffer 
energy flux deviation. Kriz and Stinchcomb [45] reported the energy flux deviation 
could be as high as 43.4 degrees for the quasilongitudinal wave and 28.7 degrees for the 
quasitransverse wave. The energy flux deviation was reported to be about 17.8 degrees 
for the pure mode transverse wave. Waves propagating in the x 2 , x 3 plane 45 degrees 
between the axes with l x = 0 have the same solutions. 

The solutions for the wave modes in an orthotropic medium will now be derived. 
The Ay matrix components are now given by 



- 134 - 


X n = c ll*f + c 66^2 + c 55^3 > 
^22 = c 66^f + c 22^2 + c 44^3' 


^33 = c 55*l + c 44*2 + c 33*3> 


^12 = ^21 = ( c 12 + c 66 ^ 1 ^ 2 » 


(B.12) 


and 


^13 = ^31 — ( c 55 + C 13)U^3> 


^23 = ^32 ~ ( c 44 + c 23)l2 l 3- 


For propagation along the x* axis, these equations reduce to 


\i — c H’ 


^22 — c 66> 


(B.13) 


^33 — c 55> 

and all other components are zero. When these are substituted into equation (2.55), the 
wave solutions can be found to be 


v, = 


_1L 

P 0 


v 2 = 


( V 

C66 ! 

P 0 


(B.14) 


and 


v 3 = 


( V— 
c 55 2 

P 0 
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with the first wave being a pure mode longitudinal wave. The second wave is a pure 
mode shear wave with particle displacements along the x 2 axis and the third wave is 

pure mode shear polarized along the X 3 axis. There is no energy flux deviation for these 
waves. 

For propagation along x 2 , the components of Xy reduce to 

\l = c 66> 

^22 = c 22> (B.15) 

= c 44> 


and all other Xy = 0 . When substituted, these yield the equation 

( C 66 ~ ^O v2 )( c 22 ~ Po y2 )( c 44 ~ Po y2 ) = 0 (B.16) 


which has the solutions 


V 1 = 


c 66 


P 0 


c 22 


li 



(B.17) 


and 


Vo = 


-44 


Po 


In this case, the first and third waves turn out to be pure mode shear waves polarized 
along the x x and x 3 axes respectively while the second wave is a pure mode compres- 
sional wave. Again, there is no energy flux deviation for these waves. 
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Now the case of wave propagation along the x 3 axis in an orthotropic material 
is presented. With 1, = 1 2 = 0, and 1 3 = 1 , the components of Xy reduce to 


“ c 55> 
^22 = c 44> 
^33 = c 33> 


(B.18) 


and all other X^ components equal to zero. These yield the solutions for the three modes 
of propagation as 


Vi = 


Vo — 


c 55 

Po 


C 44 

Po 


(B.19) 


and 


Vo = 


c 33 

Po 


In this case, the first and second waves are pure mode shear waves. The first wave is 
polarized along the x x axis. The particle displacements for the second wave are along 
the x 2 axis. The third wave is a pure mode compressional wave. There is no energy flux 
deviation for these waves. 

As with the case of transverse isotropy, the pure mode waves do not yield 
sufficient data to determine all of the orthotropic moduli. Waves must be propagated 
off axis to determine the moduli c 12 , c 13 , and c . Waves propagating in the Xj, x 2 plane 
45 degrees from each axis with 1 3 = 0 are considered first. The direction cosines of the 
wave normal are given by l x = 1 2 = , and 1 3 = 0 . Thus the Xjj matrix components 


i 
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are 


^11 = C 66 )> 

^22 — ^"( C 66 + c 22 )» 

\j3 “ "^"( c 55 + c 44). (B.20) 

^12 = ^21 = ' j ( C 12 + c 66 )> 

and all other Xy = 0 . The solutions for the waves propagating along this direction can 
then be determined to be 

P0 Y l = 7 [ C H + 2C 66 + c 22 + V( c ll ~ c 22^ + 4 (c 12 + C 66 ) 2 J 

Po y 2 = 'j( c 44 + c 55)> (B.21) 

and 

Po V 3 “ ~ J^ c 1 1 H* 2c 66 4- C 2 2 — \/( c ll “ c 22 ) 2 + 4(c 12 + 

The second wave is a pure mode shear wave while the first is quasilongitudinal and the 
third is a quasishear wave. The energy flux deviation for these waves is small because 
there is only a small deviation from isotropy in this plane. 

The next off axis direction considered is that in the x^, X 3 plane 45 degrees from 
either axis where the direction cosines for the wave normal are given by 
lj = 1 3 = and I3 = 0 . Therefore, 

^11 = + c 55)> 
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^22 ~ "^’( c 66 + c 44 )> 

^33 ~ ”^’( c 55 + c 33 )> ( B . 22 ) 

^13 = ^31 = "^"( c 13 + c 55)< 

and all other Xy = 0 . The solutions for these waves are then found by solving the equa- 
tion 

(X» - P„V 2 )[(X„ - p 0 v 2 )(X 3 3 - Po v 2 ) - Xfal - 0. (B.23) 

The three roots of this equation are given by 

Po v l ~ [ C ll + 2c 55 + c 33 + \A C 11 — C 33f + 4 ( c 13 + c 55) 2 ]> 

P0 V 2 2 = ~( c 44 + c 66)> ( B24 ) 

and 

Po v 3 2 = ~ [ c ll + 2 c S 5 + c 33 “ V ( c ll “ c 33) 2 + 4 ( c 13 + c 55)*} 

The first and third waves are quasilongitudinal and quasishear mode waves respectively 
while the second is pure mode shear. The energy flux deviations are again large as was 
the case for the similar waves for transverse isotropy. 

The final direction of propagation considered is that in the X 2 , X 3 plane 45 
degrees from either axis. This is similar to the previous case with the final solutions 
reducing to 

P0 V 1 2 = “ [ c 22 + 2c 44 + c 33 + V ( C 22 “ C 33) 2 + 4 ( C 23 + c 44)* 

Po v 2 = 7r( C 55 + c 66)> 


(B.25) 
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and 

Po V 3 - J [c 2 2 + 2c 44 + C33 - V(c 2 2 - C33) 2 + 4 (c 23 + c 44 ) 2 J 

Again, the second wave is pure mode shear while the first is quasilongitudinal and the 
third is quasishear. These equations reduce to those given for transverse isotropy if the 
appropriate substitutions are made. These equations allow the determination of all of 
the second order moduli for an orthotropic material. 
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APPENDIX C 


The relationship between the measured time between echos (t) and the true time 
of flight ((5) was given by McSkimmin [19] as 


t = p<5 — 


P7 


360F 


+ 


/ \ 

_n 

F 


(C.l) 


In this equation, p is an integer expressing the difference between the number of echos 
used for the overlap measurement. Since the first and second echos were used for these 
measurements, it was always true that p = 1 . The variable n in this expression is an 
integer which gives the number of cycles of incorrect overlap between the two echos. If 
the echos are overlapped correctly, then n = 0 . If the second echo is overlapped one 
cycle before the first then n = — 1 and if the reverse is true then n = 1 . 7 is the phase 
angle (in degrees) accounting for phase shifts at the transducer-specimen interface. F is 
the frequency of the tone burst and not the repetition rate frequency. 

Substituting the fact that p = 1 simplifies this equation to 

/ ^ ( \ 


t = 8 - 


SL 


360F 



(C.2) 




Thus both 7 and n ihust be determined to evaluate the true round trip time. Therefore 
a second equation is introduced which is 






(C.3) 


This equation expresses the change in measured time between echos as a function of a 
change in the frequency of the pulse from a higher frequency (F h ) to a lower frequency 
(Fi) . The phase shifts at these two frequencies are designated 7 h and 71 . For these 
measurements F h was chosen to be that of the resonant frequency of the transducer and 
Fj was chosen to be 0.9 times that value. Thus, this equation can be rewritten as 
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/ 






360 



(C.4) 


where F r is the resonant frequency of the transducer. 

The phase angles 7 ] and 7 r can now be calculated from standard electrical 
transmission line theory. It is defined that Z lf Z 2 , and Z 3 are the mechanical impedances 
of the transducer, bond material, and sample respectively. Z 4 is the mechanical 
impedance of the transducer-bond composite resonator and it is defined that 


r 



(0.5) 


Therefore, 


Z 4 


= iZ 2 


r tan k 2 l 2 + tan 
r — tan k^ tan k 2 l 2 ^ 


(0.6) 


where 1 ,, 1 2 are the thicknesses of the transducer and bond respectively, 
and kj, k 2 are the wave numbers for the transducer and bond given by 


k , = 


CO 

v i 


27 r 

\ 


(C.7) 


where Vj are the respective velocities of sound in the transducer and bond 
and Xj are the wavelengths of sound in the transducer and bond. 

Now it is assumed that the bond thickness is zero, and thus equation (C. 6 ) can be 
reduced to 


Z 4 = iZ x tan k 1 l 1 . 

Since at the resonant frequency of the transducer, 



( 0 . 8 ) 


(0.9) 


it is true that 



— 180 degrees. 


(CIO) 


At F! , 


Therefore at F r , 


kjlj = (0.9 • 180) degrees. 


(C.ll) 


and at F 1 , 


Z 4 = 0 


Z 4 = iZ l tan (0.9 • 180) 


(C.12) 


(C.13) 


For an elastic wave reflecting from the transducer-specimen interface, the ratio 
of the reflected pressure (P r ) to the incident pressure of the wave (P;) is given by 


Pr 

Pi 


Z 4 — Z 3 

Z 4 + z 3 


(C.14) 


and the phase angle is given by 


tan 7 = 


At F = F r , 


Pr 

Pi 


Im 

Pr 

Pj 


p r 

Re 



Pi 

^ 1 j 

0 - 

-z 3 


0 + z 3 


(C.15) 


(C.16) 


and therefore 
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P r _ ^4 ~ Z 3 

Pj iR4 + Z 3 


(C.18) 


where R 4 = tan (0.9 • 180) . 
Thus can be calculated to be 


7l = tan 1 


2Z3R4 

r 4 2 -z 3 2 ; 


(c.i8) 


Therefore 71 can be calculated from knowledge of the mechanical properties of the sam- 
ple and the transducer. 


Since 7 r = 0 , the equation for the change in measured time of flight as a func- 
tion of frequency can be reduced to 


At(0 ) = (O.Hln-7i/324) 

F r 


(C.19) 


where the superscript (0) implies that this equation is true for zero bond thickness. 
If n = 0 , this equation reduces to 




324F r 


(C.20) 


where the second superscript 0 implies that n = 0 . 

If measurements are made at F r and F { for a number of different overlap numbers 
(n) and the At's calculated for each n, the correct choice of n = 0 can be determined by 
comparing the experimental results with those predicted by the previous equation. The 
agreement will not be exact, but the one closest to the calculated value of At(°’°) and 
obeying the stipulation that 

^measured ^ ^ (C.21) 


which implies a positive bond thickness will be the correct choice for proper overlap. 
This theory can be extended to account for nonzero bond thickness. However, the 
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corrections were determined to be less than experimental error for these measurements. 
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APPENDIX D 


100 

REM 

PROGRAM 

TO COLLECT AND ANALYZE DATA FROM PU 

LSE ECHO 

118 

REM 

OVERLAP ULTRASONIC TECHNIQUE. WRITTEN BY PROSSER 

120 

REM 


9/25/85. 


150 

REM 




140 

REM 


variable list 


178 

REM 




130 

REM 

van ah] 

.e description type 


190 

REM 

f $ 

naae of file character 

\ Ji.. 1 

200 

REM 


transducer saterial 

\ 

218 

REM 

b* 

bond saterial " 

1 3 2 ) 

220 

REM 

5 $ 

saiple ID 

(32 ) 

TTfl 

i.-JL 

REM 

d$ 

date data taken ' 

(28) 

240 

REM 

ill 

data pts/over lap nuieric 


252 

REM 

m2 

series foriat of il 


25 a 

REM 

T7 

AVS. EXP. TRANSIT TIME " 



REM 

rn 

AVG. EXP. DELTA T ( 6 ) ’ 


25E 

REM 

Em, 

,F(l) TRANSIT TIME, CORR. DELTA T 

a 


REM 

F ^ 

XDUCER RES. FRE 8 . IN HZ * 


263 

REM 

ji 

I3MAX1 


261 

REM 

F3 

DATA FILE NO. 


r »in 

to- Lite 

REM 

I 

FILE LENGTH COUNTER 


pit 

fjv 

REM 

Cl 

FRES, DIVISION BY MATEC 110 

0 

244 

REM 

Pill, 

,3(1) AVS. DELTA T AT RES., LOW FREQ 

H 

245 

REM 

SCI) 

,RU> ST. DEV. OF TIME AT RES., LON 

FREQ. * 

244 

REM 

DC!) 

DELTA T S 

n 

247 

REM 

U(I) 

UNCERTAINTIES IN DELTA T 

■ 

248 

REM 

H,X 

SUM TOTALERS FOR AVG.’S 

■ 

278 

REM 

t 

saiple teip. 1 


280 

REM 

fl 

xducer resonant freq. " 


290 

REM 

1,K 

INTEGER LOOP COUNTERS ’ 


300 

REM 

1 (i ,1 

«) data at low freq. (.9fl! " 


310 

REM 

h(i,l 

■:i data at res. freq. 


320 

REM 




330 

REM 

character variable dimensioning 


340 

REM 




350 

DIM F$ 

!32) ,Xf(- 

52) ,B$ (32) ,D$(2B) ,SI (32) 


340 

REM 




370 

REM 

file naming and variables 


330 

REM 




390 

REM 





403 PRINT "*that is the naie of the file’" 

410 INPUT F$ 

420 PRINT "what is the saiple ID?’ 

430 INPUT SI 

440 PRINT " what is the xducer laterial?" 

450 INPUT Xt 

440 PRINT "(that is the bond laterial 7 * 

470 INPUT Bl 

4BB PRINT "what was the date data was taker,?’ 
490 INPUT 0* 


ORIGINAL PAGE IS 
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500 PRINT "what was the te«p. of the saeple in deg. C?* 

510 INPUT T 

520 PRINT 'what was the xducer resonant freq. in HH:?‘ 

530 INPUT Ft 

540 PRINT "what was IQHAXI?" 

550 PRINT "= aax. I of apparent echos in front of" 

560 PRINT "or behind the reference burst" 

570 PRINT "(ref. burst is taken to be the " 

590 PRINT "larger aeplitude burst)" 

590 INPUT J1 

600 PRINT "what is the * of data pts per overlap" 

610 INPUT HI 

620 PRINT "if series foreat of el was 1:2:3:4..." 

630 PRINT "then input 1" 

640 PRINT "if series foreat of e2 was 1:2:4:8..." 

650 PRINT "then input 2" 

660 INPUT H2 

670 IF H2=l THEN 710 

680 I" H2=2 THEN 710 

690 PRINT ‘expecting input of 1 or 2" 

700 6C TO 620 
710 REH 

720 REH description of echo-overlap data 

730 REH organization; each data set has 2 

740 REH arrays (1 for res. freq. and 1 for .9*res freq.) 

750 REH each with {2*!QHAXi+l>*el datapoints. 

760 REH a final data set contains xducer eaterial, 

770 REH bond eaterial, teep., etc. 

790 REH dieensioning arrays 

800 REH 

810 LET 1=2* ( Jl+i 1 -1 
820 LET K=Hl 

830 DIH H(I,K1 ,P(I) ,Q(I) ,F(I) 

840 DIH L(I,K),R(I),S(II,D(I),UU),E(n 
850 REH 

860 REH input data 

870 FOR K=1 TO HI 

880 PRINT "INPUT FREQ. DIVISION FACTOR FROH HATEC 110* 

890 INPUT Cl 

900 FOR 1=1 TO 2*(Jl+l)-l 

910 PRINT "freq. = resonant freq." 

920 Z=I-(J1+1) 

929 PRINT "PRESS RETURN TQ" 

930 PRINT "enter array aeeber";Z,K 

931 PRINT "INPUT YES TO CHAN6E FREQ. DIVISION FACTOR FROH: ";C1 

932 INPUT A$ 

933 IF AtO'YES" THEN 940 

934 PRINT "INPUT NEW FREQ. DIVISION FACTOR ’ 

935 INPUT Cl 

940 INPUT 020:H(I,K) 

941 H(I,K)=H(I,K)/10*C1 

942 PRINT H(I,K) 

943 PRINT “ " 

950 NEXT I 
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940 NEXT K 

970 FOR K=1 TO Ml 

971 PRINT "INPUT FRES. DIVISION FACTOR FROM MATED 110" 

972 INPUT Cl 

980 FOR 1=1 TO 2*(Jl+l)-l 
990 PRINT "freq. = .9«res. freq." 

1000 Z=I-(J1+1) 

1009 PRINT "PRESS RETURN TO 1 

1010 PRINT "enter array *e#her";2,K 

1011 PRINT "INPUT YES TO CHAN6E FRES. DIVISION FACTOR FROM: ";C1 

1012 INPUT At 

1013 IF AIO'YES" THEN 1020 

1014 PRINT "INPUT NEW FRES. DIVISION FACTOR* 

1015 INPUT Cl 

1020 INPUT §20:1.(1,10 

1021 L(I,K)=L(!,K)/I0 A C1 

1022 PRINT L(I,K! 

1023 PRINT ' 1 
1030 NETT I 
1040 NEXT K 

1050 REM COMPUTE AVS. 'S AND ST. DEV.'S 
1040 REM 

1070 FOR 1=1 TO 2*(Jl+l)-l 

1080 0=0 
1090 X=0 

1100 FOR K=l TO Ml 
1110 IF M2<>1 THEN 1140 
1120 0=0+K/H(I,K) 

1130 SO TO 1140 
1140 2=2 A (K-1) 

1150 0=0+Z/HiI,K) 

1140 NEXT K 
1178 P!I)=0/M1 
1180 FOR K=1 TO Ml 
1190 IF M2<>1 THEN 1220 
1200 X=X+(K/H(I,K)-P(I)! A 2 
1210 60 TO 1240 
1220 Z=2 A (K-1! 

1230 X=X+(Z/H(I,K)-P(!>! A 2 
1240 NEXT X 
1250 S(I)=(X/M1) A 0.5 
1260 NEXT I 

1270 FOR 1=1 TO 2* (Jl+1 ) -1 
1280 H=0 
1290 X=0 

1300 FOR K=1 TO Ml 
1310 IF M2<>1 THEN 1340 
1320 N=N+K/L(I,K) 

1330 SO TO 1360 
1340 Z=2 A !K-1) 

1350 0=«+Z/UI,K) 

1360 NEXT K 
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1378 8!I1=W/M1 
1380 FOR K=l TO HI 
1390 IF f!2< > 1 THEN 1420 
1400 X=X+(K/L(I ,K>-B(I) ) A 2 
1410 GO TO 1440 
1420 Z=2'!K-i! 

1438 X=X+tZ/L(I,K)-Qin)“2 
1440 NEXT K 
1458 R!I!=(X/H1!*0.5 
1460 NEXT I 
1470 REN 

1480 REN CALCULATES DELTA T AND UNCERTAINTIES IN DELTA T 
1490 REN 
1500 REN 

1510 FOR 1=1 TO 2*(Jl+l!-l 
1520 D!I)=e(I!-P!I) 

1530 U ! 1 1 =S ( 1 i +R i I) 

1548 NEXT I 
1550 REN 

1568 REN CALCULATES DELTA T(0) 

1570 REN 

1530 PRINT "NHAT IS THE INPEDANCE !*1E5 NECH. OHNS) OF” 

1590 PRINT ’THE X$; " TRANSDUCER?' 

1600 INPUT Z! 

1610 PRINT "NHAT IS THE INPEDANCE HUES NECH. OHNS) OF SAMPLE ” ; Si 
1620 INPUT Z3 

1630 PRINT "NHAT I.S THE INPEDANCE i*!E5 NECH. OHNS) OF THE ";B*;” BOND” 

1640 INPUT Z2 

1650 Z4=Zl*TANi0.9*Pi! 

1660 6=ATN i 2*Z4«Z3/ ( Z4 A 2-Z3 A 2! ! 

1661 IF 6)0 THEN 1670 

1662 G=6+PI 

1670 F2=F1*1000008 

1680 Gl=-1 / <2*P I*F2) *(6/0.9) 

1690 REN 

1700 REN IDENTIFY CORRECT N AND FIND DELTA T(0> FRON EXP. DATA 
1710 REN 
1720 REM 

1738 PRINT "COMPARE THE CALCULATED DELTA TS0! = ";G1; " (SEC.)' 

1740 PRINT "NITH THE FOLLONING VALUES" 

175B PRINT 

1760 FOR 1=1 TO 2* IJ 1+1 > -1 

1770 PRINT "DELTA T i"; I-(J1+1>; ") = ";D(D;" (SEC.)" 

1788 PRINT 
1790 NEXT I 

1800 PRINT "ENTER THE I OF DELTA T THAT MOST CLOSELY MATCHES THE" 

1810 PRINT "CALCULATED DELTA TUI" 

1820 INPUT B 
1830 N=0 
1840 REN 

1850 REN CALCULATE EXPERIMENTAL DELTA T(0) AND AVG. TRANSIT TINE 
I860 REM AND THEIR UNCERTAINTIES USING A WEIGHTED AVG. 

1870 REM 
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1880 n=0 

1390 12=0 
1900 Y1=0 
1910 Y2=0 

1920 FOR 1*1 TO 2MJ1+U-1 

1930 Z=I-(Jl+l) 

1931 N=2-B 

1940 E(I!=D(I)+N/F2*(l/0.9-l) 

1950 F ( I ) -P ( I ) +N/F2 
1960 X1=X1+E(!)M1/U(I) A 2I 
1970 Y1=Y!+F(l)*(l/S!I) A 2) 

1980 X2=X2+l/U(I) A 2 
1990 Y2=Y2+l/S(Ii A 2 
2000 NEXT I 
2010 T2=X1/X2 
2020 T3=Y1/Y2 
2030 Hl=X2 A -0.5 
2040 H2=Y2 A -0. 5 
2080 60SUB 3430 
2090 REM STORE DATA 
2100 REN 
2110 TLIST 
2120 PRINT " 

2130 PRINT "ENTER I OF (LAST) FILE TO STORE DATA ON" 
2140 INPUT F3 
2150 Z=30 

2160 2=Z+LEN (Sf ) +2+LEN (Df ! +2+LEN (B$) +2+LEN (F$) +2 

2170 Z=Z+LEN (X$) +2 

2180 Z=Z+10*(9+2»(2*!J!+1)*HU! 

2190 Z=Z+10t(6»2*(Jl+l)) 

2200 FIND F3 
2210 HARK 1,Z 
2220 FIND F3 

2230 WRITE F$,S»,Bt,DM»,Jl 

2240 WRITE T-,F1,F2,H1,H2,Z1,Z2 

2250 WRITE Z3,Z4,6,01,D,C6,C7 

2260 FOR 1=1 TO 2*(J1+1)-1 

2270 WRITE P(I),Q(I),S(I),RU),D(!),UUI 

2280 NEXT I 

2290 FOR 1=1 TO 2UJ1+11-1 
2300 FOR K=1 TO HI 
2310 WRITE H(I,K) ,L(I,K) 

2320 NEXT K 
2330 NEXT I 

2340 WRITE T2,T3,H1,H2 
2350 CLOSE 

2360 PRINT "DATA STORED” 

2370 PRINT " 

2410 PRINT 02: "FILE NAHE IS ";F$ 

2420 PRINT 02: "FILE I IS ";F3 

2430 PRINT 02: "TRANSDUCER HATER I AL IS ";X$ 

2440 PRINT 02: "BOND HATER I AL IS ",-B$ 
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2450 PRINT 82: "SAMPLE 10 IS ';Sf 

2460 PRINT 02: 'DATE DATA NAS TAKEN IS ";Df 

2461 PRINT 82: "SAMPLE TEMP NAS ";T;" C* 

2470 PRINT 82: "DATA PTS/OVERIAP IS *;M1 
2480 IF M201 THEN 2520 

2490 PRINT 02: 'SERIES FORMAT OF Ml VALUES IS;' 
2500 PRINT 02:'1:2:3:4...' 

2510 GQ TO 2540 

2520 PRINT 02: "SERIES FORMAT OF Ml VALUES IS;' 
2530 PRINT 02:'1:2:4:8...* 

2540 PRINT 02:" ’ 

2550 FOR 2=1 TO 2 

2560 IF 201 THEN 2590 

2570 PRINT 02: 'DATA AT RESONANT FREQ.* 

2580 GQ TO 2610 
2590 PRINT 02:" ' 

2600 PRINT 02: 'DATA AT . 9tRES0NANT FREQ." 

2610 PRINT 02:" ' 

2620 IMAGE 3A,1X,S 
2638 IMAGE 20,8X,S 
2640 PRINT 02: USING 2620: "Q/M" 

2650 FOR I=! TO Ml-1 
2660 IF M2<>2 THEN 2690 
2670 N=2'(I-ll 
2680 GO TO 2700 
2690 N=I 

2700 PRINT 02: USIN6 2630: N 

2710 NEXT I 

2720 IMAGE 2D 

2730 IF H2<>2 THEN 2760 

2740 H=2' (Ml-1) 

2750 GQ TO 2770 
2760 N=M1 

2770 PRINT 02: USING 2720: N 

2780 FOR 1=1 TO 2*(J1+1)-1 

2790 IMAGE 2D,1X,S 

2800 PRINT 02: USING 2790: 1-4 Jl+1 i 

2310 FOR K=1 TO Ml-1 

2820 IF ZOt THEN 2860 

2830 IMAGE 7D.2D,1X,S 

2840 PRINT 02: USING 2830: H( I, K) 

2850 GO TO 2870 

2860 PRINT 02: USING 2830: L( I, K) 

2370 NEXT K 

2880 IF ZOl THEN 2920 

2890 PRINT 02: USING 2900sH(I,Nl) 

2900 IMA6E 7D.2D 
2910 GO TO 2930 

2920 PRINT 02: USING 2900:L(I,M1) 

2930 NEXT I 
2940 NEXT Z 
2950 PRINT 02:" " 

2960 FOR Z=1 TO 2 
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2978 IF Z=2 THEN 3080 
2988 PRINT <2: 'DATA AT RESONANT FREQ.* 

2998 60 TO 3810 

3888 PRINT 82: ‘DATA AT .9+RESONANT FREQ.' 

3818 PRINT 02: ' * 

3828 PRINT 82: USIN6 3830: *8' 

3838 INA6E 2A,1X,S 

3048 PRINT 82: USIN6 3050: 'AV6. TINE (SEC.)' 

3858 INA6E 20A,1X,S 

3060 PRINT 82: USIN6 3870: 'STD. DEV. (SEC.)' 

3878 IHA6E 28A 
3088 PRINT 82:' ' 

3090 FOR 1=1 TO 2* (J 1+1 1-1 
3108 IF Z=2 THEN 3188 
3118 PRINT 82: USING 312B;!-(J1+1) 

3120 IN AGE 20,11,5 
3138 PRINT 82: USING 3148:P<!) 

3148 INA6E 10E,10X,S 
3150 PRINT 82: USING 3168:S(I) 

3168 INAGE 18E 
3178 GO TO 3210 

3180 PRINT 02: USING 3120: HJl+ll 
3190 PRINT 82: USIN6 3148:0(1) 

3200 PRINT 02: USIN6 3160:R(I) 

3210 PRINT 02:' ' 

3228 NEXT I 
3230 PRINT 02:' 1 
3240 NEXT Z 
3258 PRINT 82:' • 

3260 PRINT 02: ' IHPEDANCE OF XDUCER IS ';Z1;“ («1E5 NECH. OHHSI' 

3270 PRINT 02: "IMPEDANCE OF SAMPLE IS *;Z3;' (*IE5 NECH. OHMS)" 

3280 PRINT 82:' IHPEDANCE OF BOND IS ';Z2;' (»1E5 NECH. OHHS)' 

3290 PRINT 02: "PHASE SHIFT AT LON FREQ. IS *;6*360/<2+PI) (DEG.)' 
3300 PRINT 02: 'CALCULATED DELTA T(0) IS ';01;' (SEC.)* 

3310 FOR 1=1 TO 2»(Jl+l)-t 
3320 PRINT 02:' .' 

3330 PRINT 02: "DELTA T(";I-(J1+1);") IS ";D(I);" (SEC.)* 

3340 PRINT 02: "UNCERTAINTY CjHJl+l);') IS *;U(I) ; ' (SEC. ) ' 

3350 NEXT I 
3360 PRINT 02:' ' 

3370 PRINT 02: "CHOICE FOR CORRECT OVERLAP IS ';B 

3380 PRINT B2:’AV6. EXPERINENTAL DELTA T(N=0) IS *;T2;' (SEC.)" 

3390 PRINT 02: "ST. DEV. IS '; HI; ' (SEC.)* 

3400 PRINT 02:" • 

3410 PRINT 02: "AVG. TRANSIT TINE FRQH EXP. DATA IS *;T3;' (SEC.)" 

3420 PRINT 02: "ST. DEV. IS ';H2; ' (SEC.)' 

3421 PRINT 02:' ' 

3422 PRINT 02: "BOND ANGLE IS ";C6;' DEGREES' 

3423 PRINT 02: "TINE CORRECTION FOR BOND IS ';C7;' SEC.' 

3424 GO TO 3740 
3430 REN 

3440 REN BONO CORRECTION 
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3450 REM 

3460 DIH V(30),C(30),Q(30) 

3470 FOR 1=1 TO 30 
3480 Cl=(I-l)/2*(2*PI)/360 

3490 Z5=Z2*(Z2/Zl*TAN(Ct>+TAN(l.9*PI)!/(Z2/Zl-TAN(B.9*PI)*TAN(ClH 

3500 C2=-2*Z5*Z3 

3510 C3=Z3*Z3-Z5*Z5 

3520 C4=-2*Z2*TAN(C1)*Z3 

3530 C5=Z3*Z3-Z2*TAN (Cl 1 * ( Z2*TAN (C 1 1 ) 

3540 V( I )=ATN (C4/C5) 

3541 V(I)=V(I)*360/(2*PI) 

3550 C(I)=ATN(C2/C31 

3551 IF C(I1>0 THEN 3555 

3552 C(I)=C(I)+PI 

3555 C(I! =C(! 1*360/ (2*PI) 

3560 0 ( I ) =V ( I ) / (36B*F2)-C (I)/ (360*0. 9*F2) 

3570 NEXT I 

3580 PA6E 

3581 VIEWPORT 10,120,10,90 

3590 WINDOW I, 15,0( l )+5*Q ( 1 1 / 18 , I 

3591 HOVE 0,0 

3592 FOR 1=2 TO 16 STEP 2 

3593 HOVE 1-2,0 

3594 PRINT 'ICH’jI-Z; 

3595 NEXT I 

3600 AXIS 2, 0(11/10, e,0 
3610 HOVE 0,0(1) 

3620 FOR i=l TO 30 

3630 DRAW (1-11/2,0(11 
3640 NEXT I 

3650 HOVE 8,T2 

3660 DRAW 30, T2 

3661 HOVE 0,T2*H! 

3662 DRAW !5,T2+H1 

3663 HOVE 0.T2-H1 

3664 DRAW 15, T2-H1 

3665 HOVE 0,6*0(11/10 

3670 PRINT "INPUT PHASE A&61E OF INTERS!'' ION " 

3680 INPUT 06 

3631 HCHE 
3682 PASS 
3690 1=2*06*1 

3708 C‘=V(1! / -!3t3*F2 ( 

3710 PRINT "BOND AN3LE IS ";C6; n DEGREES * 

372? mH‘ *T>E CORRECTION FOR BOND IS “;C7;" SEC.’ 

3730 FE’l :!, i 
3748 END 
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APPENDIX E 

ORIGINAL PAGE IS 

108 REH PR06RAH TQ MEASURE STRESS VS. STRAIN 


OF POOR QUALITY 

128 REH WRITTEN BY PROSSER 7/15/86 

138 REH 

148 REH VARIABLES 

NAME 

TYPE 

158 REH 

168 REH SAHPLE ID. 

St 

CHARACTER (32) 

170 REH DIRECTION OF STRAIN l 

HI 

* (12) 

188 REH DIRECTION OF STRAIN 2 

Nt 

1 (12) 

200 REH 1 1 LOADING 

Lt 

’ (12) 

220 REH DATE 

Dt 

CHARACTER (12) 

278 REH CROSS SECTIONAL AREA OF SAHPLE (SO. IN.1A1 

« 

280 REH EXCITATION VOLTAGE FOR STRAIN GAGE l 

A2 


2?0 REH 1 ’ ' *2 

A3 


300 REH NUHBER OF OATA POINTS 

B1 

1 

330 REH V0LTA6E MEASUREMENTS FOR LOAD 

FID 

ARRAY 

340 REH VOLTAGE HEASUREHENTS FOR STRAIN 1 

V(I,1) 


350 REH ■ ■ STRAIN 2 

V(I,2) 

« 

400 REH STRAIN 1 

L(I) 

ARRAY 

410 REH STRAIN 2 

S(I) 

* 

460 REH STRESS 

D(I) 

1 

520 REH HIN STRESS 

HI 

■ 

530 REH HAX STRESS 

H2 

1 

540 REH HIN STRAIN 1 

H3 

1 

550 REH HAX STRAIN l 

04 

* 

560 REH HIN STRAIN 2 

H5 

a 

570 REH MAX STRAIN 2 

H6 

« 

571 REH SUM D(I) 

X3 

a 

572 REM SUH (D(I)) A 2 

X4 

a 

573 REH SUH L(I) 

X5 

a 

574 REH SUH (L(I)) A 2 

X6 

a 

575 REH SUH S(I) 

X7 

a 

576 REH SUH (S(I>) A 2 

X8 

a 

577 REH SUH (D(I)*l(i) 

07 

a 

578 REH SUH (D(I))*sii>) 

08 

a 

590 REH DUHHY VARIABLE 

Ct 

CHARACTER (32) 

608 REH 1 

At 

a 

630 REH DUHHY VARIABLE 

A5 

a 

648 REH INTEGER COUNTERS 

I,J 

a 

650 REH DUHHY ARRAY FOR PL0TTIN6 

P(I) 

ARRAY 

660 REH " VARIABLES FOR PL0TTIN6 

A8,A9 

NUMERIC 

670 REH PLOTTER i 

B2 

I 

680 REH REGRESSION COEF. FOR STRESS-STRAIN 1 

R1 

a 

681 REH " STRESS-STRAIN 2 

R2 


690 REH SLOPE FOR STRESS-STRAIN 1 

Cl 

a 

691 REH ' ' STRESS-STRAIN 2 

C2 

a 

692 REH Y-INTERCEPT FOR STRESS-STRAIN 1 

El 

» 

693 REH 1 ■ STRESS-STRAIN 2 

E2 

a 

700 REH 

710 REH DIMENSION CHARACTER VARIABLES 
720 REH 
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730 OIK Sf (32) ,N$(12),N$(12),L$(12) ,D$(12) 

740 DIN C$(32) ,AS(32) 

750 REM 

760 REM INPUT INITIAL VARIABLES 
770 REM 

790 PRINT 'INPUT SAMPLE ID.' 

790 INPUT SI 

800 PRINT 'INPUT DIRECTION OF LOADING" 

810 INPUT LI 

820 PRINT 'INPUT DIRECTION OF STRAIN l* 

830 INPUT NS 

840 PRINT 'INPUT DIRECTION OF STRAIN 2' 

850 INPUT NS 

860 PRINT "INPUT EXCITATION VOLTAGE FOR STRAIN GA6E 1 
880 INPUT h2 

890 PRINT 'INPUT EXCITATION VOLTAGE FOR STRAIN GAGE 2 

891 INPUT A3 

900 PRINT 'INPUT DATE" 

910 INPUT DS 

990 PRINT "INPUT CROSS SECTIONAL AREA OF SAMPLE (SB. IN.)” 

1000 INPUT A1 
1300 REM 

1310 REN DIMENSION VARIABLES MAXIMUM 200 DATA PTS. TAKEN ABOUT EVERY 
1320 REN ONE SECOND 

1330 DIN F i 200) ,V(280,2! ,L(200) ,S(200) ,D(200> 

1350 REN 

1360 REN SET UP INSTRON AND TAKE DATA 
1370 REN 

1380 PRINT 'PRESS RETURN TO BEGIN TAKING DATA" 

1390 PRINT "PRESS BREAK TO STOP TAKIN6 DATA" 

1391 PRINT 'ENTER RUN 1510 TO RESUME PR06RAH" 

1400 INPUT AS 

1410 FOR 1=1 TO 200 
1420 INPUT *24:F(I) 

1430 INPUT 812: V! 1,1) 

1440 INPUT €23: V (1,2) 

1450 B1«I 

1470 FOR J=1 TO 100 
1480 3?=J 
1490 NEXT J 
1500 NEXT I 

1510 REN CALCULATE STRESSES AND STRAINS 
1520 REN 

1540 FOR 1=1 TO Bt 

1550 D(I)=(F(I)-F!1) ) /A1 *1000*6895 

1560 L ( I ) =-4* ( V ( I , 1 ! -V ( 1 , 1 ) I / ( A2*2 . 075) 

1570 S!I)=-4*(V(I,2)-V(1,2))/(A3*2.075) 

1580 NEXT I 
1890 REN 

1900 REN CALCULATE NAX AND NIN STRESS, STRAINS 
1910 Nl«l 


-155- 


1921 M2=-l 

1930 03=1 

1940 04=-l 

1950 05=1 

I960 06*- 1 

1970 FOR 1=1 TO 01 

1930 IF L( 11)03 THEN 2000 

1990 03=1(1) 

2000 IF LUK04 THEN 2020 
2010 04=1(1) 

2020 IF D(I) >01 THEN 2040 
2030 01=D(I) 

2040 IF D (I ) <02 THEN 2060 
2050 02=D(I) 

2060 IF S ( I ) >05 THEN 2080 
2070 05=S ( I ) 

2080 IF S ( I ) <06 THEN 2100 
2090 06=5(1) 

2100 NEXT I 
2110 REN 

2120 REN LINEAR REGRESSION 

2130 REN 

2140 SOSUB 7000 

2150 REN PLOT STRESS AND STRAIN 

2160 REN 

2170 60SUB 2190 

2180 60 TO 3200 

2190 B2=32 

2200 FOR 1=1 TO 2 

2210 PAGE 

2211 PRINT 0B2:'DATE: SAMPLE: ';St;’ DIR. OF STRESS: ';L* 

2212 PRINT 882:"DIR. OF STRAIN 1: DIR. OF STRAIN 2: ';N$ 

2220 VIEWPORT 20,120,30,90 

2230 DIM P (8) 

2240 P(l)=Wl 

2250 P!2)=N2 

2260 P(3)=5 

2270 IF 1=1 THEN 2310 

2230 P(5!=N5 

2290 P(6)=W6 

2300 GO TO 2330 

2310 P(5)=N3 

2320 P(6)=W4 

2330 P(7)=5 

2340 P5=3 

2350 60SUB 5000 

2360 P5=7 

2370 60SUB 5000 

2380 WINDOW P(l) ,P(2) ,P(5) , P ! 6) 

2390 AXIS 0B2:P(3) ,P(7) ,0,0 
2400 P5=4 
2410 A$=*HHH‘ 
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2428 SQSUB £888 
2438 P5=0 
2448 Af:" 

2458 SOSUB 6880 
2468 IF 1=2 THEN 2570 
2470 HOVE 8B2:D(1),LU) 

2480 FOR J=1 TO B1 
2490 DRAW 8B2:D!J) ,L(J) 

2550 NEXT J 

2560 SO TO 2660 

2570 HOVE SB2:D(1),S(1! 

2530 FOR J=l TO B1 

2590 DRAW SB2:D(J),S(Ji 

2658 NEXT J 

2660 IF 1=2 THEN 2703 

2670 A8=«3 

2680 A9=N4 

2690 SO TO 2720 

2700 A8=N5 

2710 A9=W6 

2720 HOVE 8B2:W1,A9 

2738 PRINT 0B2: 'HHHHHJJ'; 

2740 SO TO 2810 
2750 FOR J=l TO LENiCI! 

2760 A$=SE6(C*,J,1> 

2770 PRINT *B2:A$;'HJ"; 

2788 NEXT J 
2790 HOVE 882: ill, AB 
2800 PRINT 8B2:'JJJ 
2301 SO TO 2842 
2810 IF 1=2 THEN 2840 
2820 C$='STRAIN 2’ 

2330 SO TO 2750 
2340 Cf="STRAIN l 1 
2841 SC TO 2750 
2342 Cf="STRESS (Pa) * 

2843 PRINT 8B2:C* 

2850 IF 82 0 32 THEN 2880 

2860 29=30 

2870 SO TO 2890 

2380 29=2 

23 < ’0 HQHE 

2900 FOR 28=1 TO 29 
2910 PRINT BB2:" * 

2920 NEXT 23 

2930 PRINT 882:*HIN. STRESS: ";H1; * Pa HAX. STRESS: ';1I2;' Pa’ 
2960 IF 1=2 THEN 2990 

2970 PRINT §82: 'HIN. STRAIN 1: *;W3; ‘ HAX. STRAIN 1: ’;W4 
2980 SO TO 3020 

2990 PRINT 8B2: 'HIN. STRAIN 2: , ;«5;‘ HAX. STRAIN 2: ";W6 

3020 IF 1=2 THEN 3050 

3030 PRINT §B2: 'SLOPE IS: (Pa'-ll* 
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3149 60 TO 3080 

3050 PRINT 0B2: "SLOPE IS: ";C2;" (PaMP 

3089 IF 1*2 THEN 3118 

3090 PRINT 0B2PY-INTERCEPT: ";E1 
3100 SO TO 3130 

3110 PRINT «B2: "Y- INTERCEPT: ";E2 

3130 IF 1*2 THEN 3160 

3140 PRINT ®B2: •REGRESSION CQEF. : ";Rl 

3150 60 TO 3170 

3160 PRINT 0B2: "REGRESSION COEF.: ";R2 
3170 INPUT AS 
3180 NEXT I 
3190 RETURN 

3200 HONE 

3201 PA6E 

3210 PRINT "HARD COPY OF PLOTS? YES=1 N0*2" 

3220 INPUT A5 

3230 IF A5<>1 THEN 3270 

3240 PRINT "INPUT PLOTTER NUMBER" 

3241 INPUT B2 
3250 60SUB 2200 
3260 60 TO 3300 

3270 IF A5*2 THEN 3308 

3280 PRINT "BAD CHOICE, TRY A6AIN" 

3290 60 TO 3210 

3300 PRINT 'STORE DATA? YES =1 N0=2" 

3310 INPUT A5 

3320 IF B501 THEN 3350 

3330 60SUB 3080 

3340 60 TQ 3380 

3350 IF A5=2 THEN 3380 

3360 PRINT "BAD CHOICE, TRY A6AIN" 

3370 60 TQ 3300 
3380 REN 

3460 PRINT "ANOTHER RUN? YES=l NQ=2" 

3470 INPUT A5 

3480 IF A5<>! THEN 3500 

3490 SO TQ 1370 

3500 IF A5=2 THEN 3538 

3510 PRINT “BAD CHOICE, TRY A6AIN" 

3520 60 TQ 3460 
3530 END 

5000 REN SUBROUTINE FOR PL0TTIN6 

5001 REN P(P51=NININUN 4 OF TICS 
5010 P1=(P(P5-1)-P(P5-21)/PIP5) 

5020 P2=1B A INT (L6T (PI) ) 

5030 Pl*Pl/P2 

5048 IF PI >2 THEN 5080 

5050 IF Pl=l THEN 5120 

5060 P2=2*P2 

5070 60 TQ 5120 

5080 IF PI >5 THEN 5110 
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5098 P2=5*P2 
5100 SO TO 5128 
5110 P2=10»P2 

5120 REM ADJUST DATA MINIMUM 
5130 Pl=INT(P(P5-2)/P2) 

5140 P3=P2*(Pl+2) 

5150 IF P3(P(P5-2) THEN 51B0 

5160 P3=P3-P2 

5170 SO TO 5150 

5180 P!P5-2)=P3 

5190 REM ADJUST DATA MAXIMUM 

5200 Pl=INT(P(P5-l)/P2) 

5210 P3=P2»iPl-2) 

5220 IF P<P5-1 XP3 THEN 5250 
5230 P3=P3+P2 
5240 60 TO 5220 
5250 P(P5-1I=P3 

5260 REM P(P5)=ADJUSTED TIC INTERVAL 

5270 P(P5)=P2 

5280 RETURN 

6000 REM LABLE AXIS 

6010 P4=P<P5-1) 

6020 P!4)=P(1) 

6030 PI8)=P(5) 

6040 P3-ABS(P(P5-3)+P4) MAX ABS (P (P5-2) -P4 ) 

6050 P3=INT (L6T (P3) +1 . 0E-3) 

6060 P2=10 A -P3 
6070 Pl=P (P5-2) -P4/2 
6080 P!P5)=P(P5)+P4 
6090 IF P(P5))P1 THEN 6140 
6100 MOVE 8B2:P!4) ,P!8) 

6110 PRINT At; 

6120 PRINT «B2: USING *-D.2D t S':P(P5)*P2 

6130 SO TO 6088 

6140 IF P3=0 THEN 6180 

6150 P(P5)=P1 

6160 MOVE 082:P(4),P!8) 

6170 PRINT 0B2: USIN6 •2A,+FD,S": ' E";P3 
6180 RETURN 
7008 REM 

7010 REM SUBROUTINE TO DO LINEAR REGRESSION ON DATA 

7020 REM 

7030 X3=0 

7040 X4=0 

7050 X5=0 

7068 X6=0 

7070 X7=0 

7080 X8=0 

7090 N7=0 

7100 08=8 

7110 FOR 1=1 TO B1 
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7120 X3=X3+D (I) 

7130 X4=X4+D ( I ) '2 
7140 X5=X5+L ( I ) 

7150 X6=X6+L 1 I ) A 2 
7160 X7=X7+S ( I ) 

7170 X8=X3+S(I) '2 
7180 W7=W7+D ( I ) *L ( I ) 

7190 WB=N8HMI)*S(I) 

7200 NEXT I 

7210 Ci*<Bl*4i7-X3*X5> / (B1»X4-X3“2) 

7220 C2= (B1*»8-X3*X7) / (B1*X4-X3'2J 
7230 El= ( X4*X5-X3*W7> / (B1*X4-X3 A 2) 

7240 E2= < X4*X7-X3*W8) / (B1*X4-X3 A 21 
7250 R1=(B1»W7-X3*X5>/((BIM-X3 A 2>*<S1*X6-X5 A 2)) A 0.5 
7260 R2=(BW8-X3*X7)/({BIM-X3 A 2)«<B1«X8-X7 A 2)) A 0.5 
7270 RETURN 
8000 REH 

3010 REN SUBROUTINE TO STORE BATA 

8020 PRINT "INPUT TITLE FOR STESS-STRAIN 1 FILE (MAX 64 LETTERS)* 
3030 OIN B*(64),K*(64) 

3040 INPUT Bl 

3050 PRINT ’INPUT TITLE FOR STRESS-STRAIN 2 FILE (MAX 64 LETTERS)* 
3060 INPUT Kl 
3070 A5=0 

8071 A5=A5+LEN(B$)+LEN{KI)+1B*(B1*2I+10 
8080 TLIST 
8890 PRINT *’ 

8100 PRINT ‘INPUT I OF (LAST) FILE TO STORE DATA FOR STESS STRAIN 1* 

8110 INPUT N3 

9120 FIND H3 

8130 NARK 1 , AS 

8148 FIND N3 

9150 WRITE Bl 

8160 WRITE Bl 

3178 FOR 1=1 TO Bl 

3180 WRITE D(I),L(I) 

8190 NEXT I 
8200 CLOSE 
8210 FIND N3+1 
9220 NARK 1,A5 
8230 FIND H3+1 
3240 WRITE Kl 
3250 WRITE Bl 
8260 FOR 1=1 TO Bl 
8270 WRITE D(t),S(I) 

B2B0 NEXT I 
9290 CLOSE 

3380 PRINT ‘DATA FOR STRESS-STRAIN 1 STORED ON FILE ";N3 
8310 PRINT "DATA FOR STRESS-STRAIN 2 STORED ON FILE *;N3+l 
3320 RETURN 
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APPENDIX F 


In this appendix, the derivations of the independent third order elastic moduli for 
the symmetries of transverse isotropy and orthotropy are considered. After considera- 
tions of the symmetry of the stress and strain tensors and the existence of a strain 
energy function, the number of independent third order moduli for the most general 
anisotropic solid is reduced to 56. The number of independent third order elastic 
coefficients can be reduced further if the effect of the elastic symmetry of the material is 
considered. As with the case of the linear elastic moduli, the effect of symmetry rota- 
tions on the strain energy function which is extended to include terms up to third order 
can be used to determine the independent moduli for the different symmetries. Since 
only third order terms are being considered, only the portion of the strain energy with 
third order terms (<J> 3 ) need be considered. It is given by 

$3 = tE c aaa ? ?a + ” E c aabVaVb + E c abdVa 1 7bVd- (F I) 

0 A 1 A*B A<B<D 

The process is exactly the same as was the case of reduction of the second order 
moduli. The strain energy function is written out. Then the effect of a specified rota- 
tion on the strain tensor components is determined and then substituted into the strain 
energy function. Then the rotated and original strain energy functions are equated 
term by term to determine the relationships between the moduli. However, it would be 
both too time consuming and space consuming to write out the necessary equations for 
each symmetry rotation. Thus, in this appendix only the effects on the moduli for each 
rotation will be given. These are determined by simple algebraic manipulations. 

The orthotropic moduli are determined first by considering the effects of two fold 
rotations about each of the three axes. For a 180 degree rotation about the x 2 axis 
which has a transformation matrix given by 
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-1 

0 

0 


0 0 

1 0 , 

0 -1 


(K.2) 


24 of the moduli drop out leaving 32 independent moduli. The moduli which drop out 
are c ll4 , c na , c 124 , c 12a , c 134 , c 13a , c 146 , c 16a , c^, c^, c^, c^g, c 246 , c^, 0334, 0335, c 346 , 
c 444 , c 44a , C456, c 46a , 055^, and c aaa . These can all be shown to equal zero. 

The effect of a two fold rotation about the x x axis is to reduce the number of 
independent moduli to twenty. The moduli C555, c 115 , c^, C335, c 446 , c 14a , c 126 , c 136 , 
c 236j c 340j and c 5aa all become equal to zero. This leaves the remaining independent 
moduli as c m , C 112 , c 113 , c 122 , c 123 , c 133 , c 144 , c 166; C iaa , c 222 > c 223 > c 233 > c 244 > c 256 > c 333» 

c 34 4 , C355, Csgg, and c 45a . As with the case of the second order moduli, the third sym- 
metry rotation about the x 3 axis does not reduce the number of independent moduli 
further. Thus the 20 independent moduli for an orthotropic material are those listed 
above with all other c ABD = 0 . 

Again, the transversely isotropic moduli can be determined by further reducing 
the orthotropic moduli by applying the additional symmetry condition that the material 
be isotropic in the x Xf x 2 plane. This means that the material can be rotated by any 
angle around the x 3 axis and be in an elastically equivalent position. This reduces the 
number of independent moduli to nine. They are c m , c u2 , c 113 , c 123 , c 133 , c 144 , 
c i66» c 333> an d c 344 • The following relations hold for the remaining nonzero moduli 


c 222 — c lll» c 355 — c 344 


c 223 “ c 113> c 166 “ c 266 “ ^ ( c lll c 112) 


c 255 ““ c 144> c 223 “ c 133 ( F -3) 

c 244 ~ c 155> c 366 = ~( c 113 ~ c 12s) 
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C 122 “ c 1 12* c 456 “ 9 ( c 155 c 144)- 
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APPENDIX G 

100 REM PRQ6RAM TO MEASURE P2L2 FREQUENCY 

110 REM CHANGES AS A FUNCTION OF APPLIED LOAD AND PLOT AND 

120 REM STORE DATA. WRITTEN BY PROSSER 2/25/86 


130 REM 

140 REM VARIABLES 

NAME 

TYPE 

150 REM 

160 REM SAMPLE ID. 

SI 

CHARACTER (32) 

170 REM MODE OF NAVE 

Ml 

“ (12) 

180 REM DIRECTION OF PROPAGATION 

Nl 

* (12) 

190 REM • ' ' POLARIZATION 

PI 

' (12) 

200 REM ' * LOAD I NS 

LI 

' (12) 

220 REM DATE 

Dl 

CHARACTER (12! 

270 REM CROSS SECTIONAL AREA OF SAMPLE (SO. IN. 

)A1 

■ 

300 REM NUMBER OF DATA POINTS 

B1 

■ 

330 REM FREQUENCY MEASUREMENTS (Hi) 

F (I ) 

ARRAY 

340 REM VOLTAGE MEASUREMENTS FOR LOAD (VOLTS) 

VII) 

■ 

400 REM STRESS (Pa) 

L ( I ) 

ARRAY 

460 REM DELTA F/F AT EACH POINT 

0(1) 

■ 

520 REM MINIMUM DELTA F/F 

Nl 

■ 

530 REM MAXIMUM DELTA F/F 

N2 

> 

540 REM MINIMUM STRESS (Pa) 

N3 

> 

550 REM MAXIMUM STRESS 

N4 

1 

571 REM SUM D ( I ) 

X3 

I 

572 REM SUM !D(I)) A 2 

X 4 

■ 

573 REM SUM L(I! 

X5 

■ 

574 REM SUM (L(I)P2 

X6 

■ 

577 REM SUM tO(I)*l(i! 

N7 

1 

590 REM DUMMY VARIABLE 

Cl 

CHARACTER (32 

600 REM ' 

At 

• 

630 REM DUMMY VARIABLE 

AS 

0 

640 REM INTE6ER COUNTERS 

M 

l 

650 REM DUMMY ARRAY FOR PLOTTING 

P(I) 

ARRAY 

660 REM " VARIABLES FOR PLOTTING 

AB,A9 

NUMERIC 

670 REM PLOTTER 1 

B2 

l 

600 REM RE6RESSI0N COEF. FOR STRESS-DELTA F/F 

R1 

* 

690 REM SLOPE FOR STRESS-DELTA F/F 

Cl 

■ 

692 REM Y-INTERCEPT FOR STRESS-DELTA F/F 

El 

a 


700 REM 

710 REM DIMENSION CHARACTER VARIABLES 
720 REM 

730 DIM St (32) ,NI(12>,Nt(12),Pt(12) ,LH12) ,DI(12! 
740 DIM C*(32),At(32) 

750 REM 

760 REM INPUT INITIAL VARIABLES 
770 REM 

790 PRINT 'INPUT SAMPLE ID.' 

790 INPUT St 

800 PRINT 'INPUT MODE OF NAVE' 

810 INPUT Ml 

B20 PRINT 'INPUT DIRECTION OF POLARIZATION* 
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S3B INPUT PI 

840 PRINT ‘INPUT DIRECTION OF PROPAGATION" 

850 INPUT N$ 

860 PRINT ‘INPUT DIRECTION OF LOADING’ 

870 INPUT LI 

900 PRINT ‘INPUT DATE" 

910 INPUT D4 

990 PRINT ‘INPUT CROSS SECTIONAL AREA OF SANPLE (SB. IN.)" 

1000 INPUT A1 
1300 REN 

1310 REN DIMENSION VARIABLES HA XI HUH 200 DATA PTS. TAKEN ABOUT EVERY 

1320 REH ONE SECOND 

1330 DIH F (200) , V (200) ,L!200) ,D(200) 

1350 REH 

1360 REH SET UP INSTRON AND TAKE DATA 
1370 REH 

1380 PRINT ‘PRESS RETURN TO START TAKING DATA" 

1390 PRINT "PRESS BREAK TO QUIT TAKIN6 DATA" 

1391 PRINT ‘THEN ENTER RUN 1510* 

1400 INPUT A$ 

1410 FOR 1=1 TO 200 
1420 INPUT «3:F(I! 

1430 INPUT 024: VII) 

1450 81=1 

1470 FOR J=1 TO 85 
1480 89=J 
1490 NEXT J 
1500 NEXT I 

1510 REH CALCULATE STRESSES AND DELTA F/F 

1520 REH 

1540 FOR 1=1 TO B1 

1550 L(I)=(V(I)-V(1))/A1«6B95*1000*-1 
1570 D(I)=(F(I)-F(1))/F(1) 

1580 NEXT I 
1390 REH 

1900 REH CALCULATE HAX AND HIN STRESS AND DELTA F/F 

1910 Nl=l 

1920 N2=-l 

1930 N3=l 

1940 N4=-l 

1970 FOR 1=1 TO B1 

1980 IF L ( I ) >N3 THEN 2000 

1990 M3=L(I) 

2000 IF L(IKH4 THEN 2020 
2010 N4=L(I) 

2020 IF DCION1 THEN 2040 
2030 U1=D(I) 

2040 IF D(IXN2 THEN 2100 
2050 N2=D(I1 
2100 NEXT I 
2110 REH 

2120 REH LINEAR REGRESSION 
2130 REH 
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2140 SOSUB 7000 

2150 RES PLOT STRESS VS. DELTA F/F 

2160 RES 

2170 SOSUB 2190 

21B0 60 TO 3200 

2190 B2=32 

2210 PA6E 

2211 PRINT 082: 'DATE: SASPLE: ’iSt;’ NAVE NODE: ';S$ 

2212 PRINT SB2: "DIR. OF PROPAGATION: ";Nt;" DIR. QF STRESS: *;L« 

2213 PRINT 882: 'DIR. OF POLARIZATION: ‘;P$ 

2220 VIEWPORT 20,120,30,90 

2230 DIS P (8) 

2240 P(5)=W1 
2250 P(6) s W2 
2260 P (3) =5 
2310 P(1)=W3 
2320 P!2)=W4 
2330 P (7> =5 

2340 P5=3 

2350 SOSUB 5000 
2360 P5 S 7 
2370 SOSUB 5000 

2380 WINDOW P < 1 ) , P (2) ,P(5I ,P(6) 

2390 AXIS 0B2:P(3) , P (7) ,0,0 
2400 P5M 
2410 AI="HHH* 

2420 SOSUB 6000 
2430 P5=8 
2440 A*»" 

2450 SOSUB 6000 
2470 SOVE §B2:L(1),D(1) 

2480 FOR J=t TO Bl 
2490 ORAN 0B2:L(J> ,D(J) 

2550 NEXT J 

2670 A8=W3 

2680 A9=N4 

2720 SOVE §B2:A8,N2 

2730 PRINT 0B2: "HHHHHJJ'j 

2740 C$ S 'DELTA F/F"’ 

2750 FOR J*1 TO LEN(CI) 

2760 AI=SES!CI,J,1) 

2770 PRINT 0B2:At; ’HJ'; 

2780 NEXT J 
2790 SOVE SB2:A8,W1 
2800 PRINT 0B2:'JJJ 
2820 Cf='STRESS lPai " 

2341 PRINT 0B2:C$ 

2850 IF B2<>32 THEN 2880 
2860 Z9=30 
2870 SO TO 2890 
2880 Z9*2 
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2890 HOME 

2900 FOR 28=1 TO Z9 
2910 PRINT 082;’ " 

2920 NEXT Z8 

2930 PRINT 082: "MIN. DELTA F/F: Wl;" MAX . DELTA F/F: ';«2 

2970 PRINT 0B2:”MIN. STRESS: ";W3; 1 (Pal";' MAX. STRESS: ”;N4;" !Pal" 

3030 PRINT 0B2: 'SLOPE IS: "jCl;" <Pa A -l>" 

3090 PRINT €B2: a Y- INTERCEPT: “;E1; ’ (Pa'-l)' 

3140 PRINT 082: "REGRESSION COEF. : ";R1 
3170 INPUT Al 
3190 RETURN 

3200 HONE 

3201 PAGE 

3210 PRINT "HARD COPY OF PLOTS? YES=1 N0=2* 

3220 INPUT A5 

3230 IF A5<)1 THEN 3270 

3240 PRINT "INPUT PLOTTER NUMBER* 

3241 INPUT 82 
3250 6GSUB 2210 
3268 GO TO 3300 

3270 IF A5=2 THEN 3300 

3280 PRINT "BAD CHOICE, TRY AGAIN" 

3290 GO TO 3210 

3300 PRINT "STORE DATA? YES =1 NQ=2’ 

3310 INPUT A5 

3320 IF A5<>! THEN 3350 

3330 SQSUB 3000 

3340 GO TO 3330 

3350 IF A5=2 THEN 3330 

3368 PRINT "BAD CHOICE, TRY AGAIN" 

3370 GO TO 3300 
3380 REN 

3460 PRINT "ANOTHER RUN? YES=1 NO=2" 

3470 INPUT A5 

3480 IF A5<>1 THEN 3500 

3490 SO TO 1370 

3500 IF A5=2 THEN 3530 

3510 PRINT "BAD CHOICE, TRY AGAIN” 

3520 GO TO 3460 
3530 END 

5000 REN SUBROUTINE FOR PLOTTING 

5001 REN P!P5) MINIMUM I OF TICS 
5010 ?t=fP(P5-l)-P(P5-2>)/P!P5) 

5020 P2=10'*INT!L8T (PI) ) 

5030 Pt=Pi/P2 

5040 IF PI >2 THEN 5080 

5050 IF Pl=t THEN 5120 

5060 P2=2*P2 

5070 60 TO 5120 

5080 IF PI >5 THEN 5110 

5090 P2=5*P2 
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5100 SC TO 5120 
5110 P2=10*P2 

5120 REN ADJUST DATA NININUH 
5130 P1=INT (P(P5-2)/P2) 

5140 P3=P2*(Pi+2) 

5150 IF P3<P(P5-2) THEN 51S0 

5160 P3=P3-P2 

5170 SO TO 5150 

5180 P!P5-2)=P3 

5190 REN ADJUST DATA NAXINUN 

5200 P1=INT !P!P5-1)/P2) 

5210 P3=P2»(Pl-2) 

5220 IF P(P5-1KP3 THEN 5250 
5230 P3=P3+P2 
5240 SO TO 5220 
5250 P(P5-1!=P3 

5260 REN P <P5) -ADJUSTED TIC INTERVAL 

5270 P ( P5> =P2 

5280 RETURN 

6000 REN LABLE AXIS 

6010 P4=P;P5-1) 

6020 P(4)=P(1) 

6030 P(8)-P(51 

6040 P3=ABS(P(P5-31+P4! NAX ABS(P'P5-2)-P4) 

6050 P3=INT!LST(P3)+1.0E-8) 

6060 P2=10"'-P3 
6070 P1=P <P5-2) -P4/2 
6080 PIP5)=PIP5)+P4 
6090 IF P!P5) >P1 THEN 6140 
6100 NOVE §B2:P(4) ,P (8) 

6110 PRINT A*; 

6120 PRINT *82: USING *-D.2D,S':P«P5}«P2 

6130 SO TO 6080 

6140 IF P3=B THEN 6180 

6150 P (P5) =P1 

6160 HOVE 0B2:P(4) ,P(8) 

6170 PRINT 0B2: USIN6 "2A,+FD,S": " E";P3 
6180 RETURN 
7000 REN 

7010 REN SUBROUTINE TO DO LINEAR REGRESSION ON DATA 

7028 REN 

7030 X3=0 

7040 1(4=0 

7050 X5=0 

7068 X6=8 

7090 N7=0 

7110 FOR 1=1 TO B1 
7120 X3=X3+D!I) 

7130 X4=X4+D(I)'2 
7140 X5=X5+L(I) 

7150 X6=X6+L(I) A 2 
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7180 B7«H7*DU)*L(I1 
7200 NEKT I 

7210 Cl= (Bl*ti7-X3*X5) / (B1*X6-X5 A 2) 

7230 E1=(X6*X3-X5#»7)/<B1*X6-X5 A 2) 

7250 Rl« (B1»M7-X3*X5) / ( (B1*X4-1(3 A 2) « (Bi«6-X5 A 2) ) A 0.5 
7270 RETURN 

8000 REN 
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APPENDIX H 

100 REN NONLINEAR DATA ANALYSIS PR06RAH (TRANSVERSE ISOTROPY) 
110 REN DETERNINE 7 OF 9 NONLINEAR COEF. BY LEAST SQUARES FIT 


120 REN USIN6 24 SAC HEASUREHENTS 
130 REN WRITTEN BY PROSSER 11-10-86 
140 REN 

150 REN VARIABLES 

160 REN LINEAR ELASTIC STIFFNESSES K(I) 

170 REN ' ■ CQNPLIANCES S(I) 

180 REN SAC NEASURENENTS H(I,1) 

190 REN UNCERTAINTIES IN SAC'S H(I,2! 

200 REN NONLINEAR COEF'S C(I,1) 

210 REN NAX. UNCERTAINTY IN NONLINEAR COEF’S C(I,2) 

220 REN PROB. * . • " C(I,3) 

230 REN RECALCULATED SAC’S H(I, 3) 

240 REN NAX. UNCERTAINTY IN RECALCULATED SAC’S HU, 4) 

250 REN PROB. 1 " HU, 5! 

260 REN CQEFICIENTS OF EQ’S (24 Efl., 7 UNKNOWNS) A(I,J) 

270 REN AU,J)*AU,N> E(I,Ji 

280 REN INVERSE OF E(I,J) FU,J) 

290 REN PRODUCT OF E AND F F1!I,J) 

300 REN PRODUCT OF F AND E F2!I,J) 

310 REN CONDITION NUNBER OF E Z3 

320 REN H( 1,1) -•’CONSTANT TERNS YU) 

330 REN INTE6ER COUNTERS I,J,L,N 

340 REN PERCENT DIFFERENCE BETWEEN NEASUREE AND CALC SAC'S VU) 
350 REN ” UNCERTAINTY IN CALCULATED SAC’S T(I,2> 

1000 REN INPUT DATA 

1001 DIN K(3) ,S<4) ,H(24,5) ,A(24,7! ,C(7,3) 

1010 PAGE 


1020 PRINT "INPUT LINEAR ELASTIC STIFFNESS HODULI (6Pa)" 

1030 PRINT "11-1,12-2,44-3" 

1040 FOR 1=1 TO 3 
1050 PRINT "INPUT C"; I 
1060 INPUT K(I) 

1070 NEXT I 
1080 PA6E 

1090 PRINT "INPUT LINEAR ELASTIC CONPLIANCE NODULI (6PA) A -1* 
1100 PRINT "11-1,12-2,13-3,33-4" 

1110 FOR 1=1 TO 4 
1120 PRINT "INPUT S";I 
1130 INPUT 5(1) 

1140 NEXT I 
1150 PA6E 

1160 PRINT "INPUT SAC’S AND UNCERTAINTIES (GPa)'-l" 

1170 PRINT "122-1,211-2,322-3,311-4,121-5,212-6,312-7,321-8" 

1180 FOR 1=1 TO 8 

1190 PRINT "INPUT SAC "; I 

1200 INPUT H(I ,11 

1210 PRINT "INPUT UNCERTAINTY" 

1220 INPUT HU, 2) 

1230 NEXT I 


ORIGINAL PAGE IS 
OF POOR QUALITY 
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124B PAGE 

1250 PRINT “SAC'S CQNT'D.' 

1260 PRINT "123-9, 213-10, 323-11, 313-12,111-13, 232-14, 132-15,231-16" 

1270 FOR 1=9 TO 16 

1280 PRINT "INPUT SAC ";I 

1290 INPUT H(I, 1) 

1300 PRINT "INPUT UNCERTAINTY" 

1310 INPUT H ! I ,2) 

1320 NEXT I 
1330 PAGE 

1340 PRINT “SAC’S CDN'T," 

1350 PRINT "011-17,022-18,012-19,021-20,013-21,023-22,031-23,032-24" 
1360 FOR 1=17 TQ 24 
1370 PRINT "INPUT SAC ";I 
1380 INPUT H ( I , 1 ) 

1390 PRINT "INPUT UNCERTAINTY" 

1400 INPUT HI I, 21 
1410 NEXT 1 

1420 REN CALCULATE COEF'S OF EQ.S 
1430 A=0 

1440 A(l,l)=S(2)/(2*K(ll) 

1450 A(1,2)=S(1) / (2»K(1) 1 
1460 A(1,3)=S(31/(2»K(1M 
1470 A(3,1)=A(1,3) 

1480 A(3,2)=A(1,3) 

1490 A(3,3)=SI4)/(2*Kil)> 

1500 A(5,l)=(S(l)+S(2))t0.25/(K(l)-K(2)l 
1510 A(5,2)=-A(5,l) 

1520 A!5,3)=0.5*S(3)/(K(1)-K(2)1 

1530 A(5 , 4) =-A (5,3) 

1531 A(7,1)=A(5,3) 

1532 A(7,2i=A(5,4i 

1540 A(7,3)=0.5*S(4)/(K(1)-K(2)1 
1550 A(7,4i=-A(7,3) 

1560 A(9,5>=S(11/(2*K(3)1 
1570 A(9,6)=S(2)/(2*K(3) ) 

1580 A!9,7)=S(3)/(2*K(31) 

1590 A(11,5)=A(9,7) 

1600 A(11,61=A(9,7) 

1610 A(ll,7!=S(4)/(2*K(3)) 

1620 A(13,5)=A(9,6) 

1630 A(13,6)=A(9,5) 

1640 A(13,71=A(9,7) 

1650 A(15,5)=A(9,5) 

1660 A(15,6)=A(9,6) 

1670 A(15,7)=A(9,7) 

1680 B=Sil)+S(2)+S(3) 

1690 D=2*S(3)+S(4) 

1700 A(17,11=B/ I2*K(1I) 

1710 A(17 f 2)=AU7,l) 

1720 A(17,31=D/(2*K(1) 1 
1730 ft 1 19, 1 ) =0. 5*B/ (KU1-KI2) ) 

1740 A (19 , 2) =-A (19,1) 
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1750 A(19,3)=D/(K(1)*K(2) 1*0.5 
1760 A ( 19,4) =— A (19,3) 

1770 A(21,5)=B/(2*K(31) 

1780 A (21 , 6) =A (21 , 5) 

1790 A(21,7)=D/(2*K(3) ) 

1800 A (23,51 =A (21 ,5) 

1810 A (23,6 ) =A (21 ,51 
1820 A (23,71 *A (21 ,71 
1830 FOR J=1 TO 12 
1840 FOR 1=1 TO 7 
1850 A(2*J,I1»A(2*J-1,I1 
1860 NEXT I 
1870 NEXT J 

1880 REM ADD CONSTANT TERNS TO SAC'S 

1890 DIN Y (241 

1900 Y(11=-(HC1,11+S(21) 

1910 Y!2)=-(H!2,1)+S(2)I 
1920 Y(3I=-(H(3, 11+8(31) 

1930 Y(41*-(H<4, 11+S(3) 1 
1940 Y (51=-(H(5, 1 1 +S( 1 ) 1 
1950 Y(6)=-(H(6,11+S(l) 1 
1960 Y(7)=-(H(7,11+S(3>) 

1970 Y(8)=-(H(8,1)+S(3) 1 
1980 Y(9)=-(H(9,1)+S(3)1 
1990 Y(10!=-(H(10,1)+S(31) 

2000 Y(ll)»-(H(il, 11+8(4)1 
2010 Y(12)=-(H(12,1)+S(41 1 
2028 Y(13)>-(H(13,1)+S(1)) 

2030 Y(i41*-(H(14,1)+S(ll) 

2040 Y(15)*-(H(15, 11+8(21) 

2050 Y(16)=-(H(16 I 1)+S(2)1 
2060 Y(17)=-(H(17,l)+8+l/(2*K(l)l) 
2078 Y(18I=-(H(18,1)+B+1/(2*K<1) 1 1 
2080 Y(191 a -(H(l9,l)+B+l/(K(l)-K(2)l) 
2090 Y (20)»-(H(20, 1 ) +B+1/ (K (1 1 -K (2) 1 1 
2100 Y(21) S -(H(21 I 1)+D+1/(2*K(31»1 
2110 Y(22)=-lH(22,l)+D+l/(2*K(3) ) 1 
2120 Y(23)*-(H(23,1)+B+1/(2»K(3))I 
2130 Y(24)=-(H(24,ll+B+l/(2+K(3)>) 
2140 REN CONFUTE A ( I , N) #A( I , J 1 
2150 DIN E (7 , 7) 

2160 FOR N=1 TO 7 

2170 FOR J=1 TO 7 

2180 E(N,J)=0 

2190 FOR 1=1 TO 24 

2200 E(H,J)=E(N,J)+A(I,H1*A(I,J) 

2210 NEKT I 
2220 NEXT J 
2230 NEXT N 

2240 REN COHPUTE INVERSE OF E(H,J) 
2250 DIN F(7,7) ,F1(7,7) ,F2(7,7) 

2260 F=INV(E) 

2270 REN CONPUTE CONDITION NUNBER 
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2280 Z1=-1.0E+100 

2290 Z2=-1.0E+100 

2300 FOR 1=1 TO 7 

2310 FOR J=1 TO 7 

2320 IF Ed.JXZl THEN 2340 

2330 Z1=E(I,J) 

2340 IF F ( I , J ) < Z 2 THEN 2360 

2350 Z2=F(I,JI 

2360 NEXT J 

2370 NEXT I 

2380 Z3=A8S(Z2*Zl) 

2390 F1=E MPY F 
2400 F2=F NPY E 

2440 REM COMPUTE NQNLINEftR COEF'S 

2450 FOR 1=1 TO 7 

2460 C(!,l)=0 

2470 C!I,2!=0 

2480 C(I,3)=0 

2490 FOR J=l TO 24 

2500 FOR L=1 TO 7 

2510 C(I,l)=C(I,l)+F(L,Il*A(J,L)*Y(J) 

2520 C(I,2)=C(I,2I+ABS(F(L,I)«A(J,L))«H<J,2I 
2530 NEXT L 
2540 NEXT J 

2550 C(!,3)=C(I,2)/23 A 0.5 
2560 NEXT I 

2570 REM RECALCULATE SAC VALUES USIN6 DERIVED VALUES 

2580 H(l,3)=-(2*K(l)*S<2)+S(l)*C(2,l)+S!2!*C(l,t)+S(3)*Cl3,i)l 

2590 H(1,3)=H(1,3)/(2*KU1) 

2600 H(3 t 3>=-(2*K(n*S(3)+S(3)*(C(l , 1)+C(2,1) ) +S ( 4 ) *C (3, 1) ) 

2610 H(3,3)=H(3,3I/(2«K(1)1 

2620 H(5,3)=(K(l)-K(2))*SU)+0.25*(S(U+S(2))*(C(i,l)-C(2,i)) 

2630 H(5,3)=H(5,3)+0.5*S(3i*(C(3,l>-C(4,l)> 

2640 H(5,3)=-H(5,3)/(K(l)-K(2!> 

2650 H(7,3)=(Kll)-K(2) >*S(3)+0.5*S(3)*(C(1,1)-C(2,1!) 

2660 H(7,3)=H(7,3)+0.5*S(4)*(C(3,t)-C(4,l)) 

2670 H (7 ,3>=-H (7,3>/(K(l)-K(2>) 

2680 H(9,3)=-(2*K(3)*SI3)+SU)*C(5,l>+S(2i*C(6,l)+S(3)*C(7,l))/(2*K(3l> 
2690 H(tl,3)=-(2*K(3)«S(4)+S(3)*<C(5,l)+C(6,l))+S(3)*C(7,l))/(2«K(3)) 
2700 H(13,3)=-(2*K (31*5(11 +S(l)*C(6,l)+S(2)*C(5,l)+S(3)*C (7,1) )/(2*K(31) 
2710 H(l5,3)=-(2*K(3)*S(2)+S(l)*C!5,ll+S!21*C(6,l)+S(3)*C(7,ll)/(2*K(3)l 
2720 H(17,3)=-(l+2»K(t)«B+B*(C(l,l)+C(2,l))+D*C(3,l))/(2*K(i)l 
2730 H(l9,3)=l+(K(l)-K(2)l*B 

2740 H(19,31=H(19,3)+0.5*(B*(C(1,1)-C(2,1H+D*(C(3,1)-C(4,1))1 
2750 H(19,3)=-H(19,31 /(K(l)-K(2) 1 

2760 H(21,3)=-(1*2*K(3)*D+B*(C(5,1)+C(6,11 )+D*C(7,ll )/(2*K(3) 1 
2770 H(23,3)=-(l+2*K(3)tB+B*(C(5,!l+C(6,l))+D*C(7,li 1 /(2*K(3) 1 
2780 FOR 1=1 TO 12 
2790 H(2*I,3)=H(2*I-t,3> 

2800 NEXT I 
2810 DIM H (7,2) 

2820 FOR 1=1 TO 2 
2830 FOR J=1 TO 7 
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2640 

H ( J, I >=ABS (C ( J 

,I+i)/C(J,l>) 

2850 

NEXT J 


2860 

NEXT I 


2870 

FOR 1=1 TO 2 


2880 

FOR J*1 TO 4 


2890 

H(J,I+3)=(0<1, 

11+0(2, 11+0(3, I))*H(J, 3) 

2900 

NEXT J 


2910 

FOR J=5 TO 8 


2920 

H(J,I+3)=(N(l, 

I)+0(2,I)+0(3,I)+0(4,I1)*H(J,3) 

2930 

NEXT J 


2940 

FOR J=9 TO 16 


2950 

H(J,T+3)=(0(5, 

l)+0(6, 11+0(7, DMHU, 3) 

2960 

NEXT J 


2970 

FOR J=17 TO 18 


2980 

H(J,I+3)=(H<1, 

11+0(2, 11+0(3, IH»H(J, 31 

2990 

NEXT J 


3000 

FOR J=19 TO 20 


3010 

H(J,I+3)=(N(1, 

11+0(2,11+0(3, 11+0(4, I1)*H(J, 3) 

3020 

NEXT J 


3030 

FOR J=21 TO 24 


3040 

H(J,I+3>=(0(5, 

11+0(6, 11+0(7, I))*H(J, 3) 

3050 

NEXT J 


3060 

NEXT I 


3070 

FOR 1=4 TO 5 


3080 

FOR J=1 TO 24 



3090 H(J,I)=ABS(H(J,I! 1 
3100 NEXT J 
3110 NEXT I 
3120 DIN V (24) 

3130 FOR 1=1 TO 24 

3140 V ( I ) =ABS ( !H ! I , 1 ) -H ( I , 3) ) /H ( I , 1 ) 1 
3150 V(I)=V(I)*100 
3160 NEXT I 
3170 DIM T (24,2) 

3180 FOR 1=1 TO 2 
3190 FOR J=l TO 24 

3200 T(J,I)=H(J, 1+3) /H(J, 31*100 

3201 T(J,I)=ABS(T(J,D) 

3210 NEXT J 

3220 NEXT I 

3230 REM PRINT OUT DATA 
3240 REM 
3250 N=32 
3260 PAGE 
3270 GQSUB 3380 

3280 PRINT "PRINTED COPY? YES =1 NO =2 " 

3290 INPUT 75 

3300 IF Z5=l THEN 3340 

3310 IF Z5=2 THEN 4010 

3320 PRINT "BAD CHOICE" 

3330 60 TO 3280 

3340 PRINT "INPUT PRINTER »* 

3350 INPUT N 
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J360 60SUB 3380 
3370 6Q TO 4010 

3380 PRINT 0N: "NONLINEAR DATA ANALYSIS FOR TRANSVERSE ISOTROPY MODEL" 
3390 PRINT «N:" " 

3400 PRINT 0N:* INPUT DATA" 

3410 PRINT IN:" * 

3420 PRINT §N: "LINEAR ELASTIC MODULI" 

3430 PRINT SN:” " 

3440 PRINT 0N:’C11 = ";K(D;" (6Pa> Sll = ";S(1>;" (GPa) A -l" 

3450 PRINT 0N:"C12 * ";K(2) ; " (6Pa» S12 = ";S!2);" (SPa)'-l" 

3460 PRINT §N: "C44 = ";K(3>;" (GPa) S13 = ";S(3);" (GPa)M’ 

3470 PRINT 0N: ’ S33 = ";S(4);" (GPa)'-l" 

3480 PRINT 0N: * " 

3490 PRINT 0N: 'VELOCITY DERIVATIVE DATA 16Pa) A -l" 

3500 PRINT «N: " * 

3510 PRINT §N:’l. H122 = ";H(1,1);’ +/- ";HU,2) 

3520 PRINT 0N: "2. H211 = ";H<2,1>;" +/- ";H(2,2) 

3530 PRINT 0N:"3. H322 = ";H(3,1>;" +/- ";H(3,2) 

3540 PRINT «N:"4. H311 = ’;H(4,1);" +/- ";H(4,2) 

3550 PRINT 0N:"5. H121 = ";H!5,D;" +/- ";H(5,2) 

3560 PRINT «N:"6. H212 = ";H(6,11;" +/- ";H(6,2) 

3570 PRINT 0N:"7. H312 = ";H(7,1!;" +/- ";H<7,2) 

3580 PRINT §N: "8. H321 = ";H(8,D ; " +/- ";H(8,2) 

3590 PRINT 0N: "9. H123 = ";H(9,1);" +/- ";H!9,2> 

3600 PRINT 0N:"10. H213 = ’;H(10,l);" ♦/- *;H(10,21 

3610 PRINT ®N: "11. H323 = ■jHdl.l);" +/- ";H(ll,2) 

3620 PRINT 0N: "12. H313 = ";H(12,l!j" +/- *;H(!2,2) 

3630 PRINT 0N:"13. H131 = ";H!13,1);," +/- ";H(13 t 2) 

3640 PRINT «N:"14. H232 = ";H(14,D; " +/- ";H(14,2) 

3650 PRINT «N:"15. H132 = ";H(15,D;" +/- ";H(15,2) 

3660 PRINT ®N: "16. H231 = "^(16,1);" +/- ";H(16,2) 

3670 PRINT 0N; *17. H011 * ";H(17,1);" +/- ";HU7,2) 

3680 PRINT «N:"18. H022 = ";H>:18, 1) ; " +/- ";H(18,2! 

3690 PRINT 0N:"19. H012 = " ; H ( 19 , 1 > ; * +/- ";H(19,2» 

3700 PRINT «N:"20. H021 = "iH(20,l);" +/- "jH(20,2! 

3710 PRINT 0N:"21. H013 = ";Hi21,D;" +/- ";H(2t,2} 

3720 PRINT 0Ni"22. H023 = ";H(22,D;" +/- ”;H(22,2) 

3730 PRINT §N:"23. H031 = ";H(23,D;* +/- ";H(23,2) 

3740 PRINT SN: "24. H032 = ";H(24 t l>;" ♦/- ";H(24,2i 

3750 PRINT 0N:" " 

3760 PRINT 0N: "THIRD ORDER ELASTIC CONSTANTS DERIVED FROM LEAST" 

3770 PRINT 0N: "SQUARES FIT AND MAX. ERROR AND PROB. ERROR (GPa)* 

3780 PRINT 0N:" * 

3790 PRINT 0N:"C111 = ";C(1,1);" +/- ";C!1,2);" +/- ";C(1,3) 

3800 PRINT 0N: "Cl 12 = ";C(2,1I;" +/- ";C(2,2);" +/- “;C(2,31 

3810 PRINT ®N: "C 1 13 = *;C(3,1); D +/- ";C(3,2>;" +/- ";C(3,3) 

3820 PRINT SN: "C123 = ";C(4,D;" +/- ";C(4,2);" +/- ";C(4,3I 

3830 PRINT 0N:"C144 = ";C(5,1I;" +/- *;C(5,2);" +/- ";C(5,3) 

3840 PRINT §N; "C155 = ";C(6,1);" +/- ";C(6,2);" +/- “;C(6,3) 

3850 PRINT 0N:"C344 = ";C(7,D;" +/- •;C(7,2!;" +/- ";C(7,3) 

3860 PRINT 0N:"C133,C333 UNDETERMINED DUE TO LACK OF DATA" 

3870 PRINT QN:’ " 

3880 PRINT §N: "RECALCULATED VELOCITY DERIVATIVES USIN6 FITTED TOEC'S" 
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3898 PRINT «N:'AND NAT. AND PRQB. ERRORS IN SANE ORDER AS BEFORE' 
3980 PRINT §N: " ' 

3910 FOR 1=1 TO 24 

3920 PRINT 0N:I;*. ';H(I,3);' ♦/- , 5H(1,4)j* ♦/- ';H(I,5> 

3930 NETT I 
3948 PRINT 0N:' ' 

3950 PRINT §N: "PERCENT DIFFERENCES, NAT. AND PROB. UNCERTAINTIES' 
3960 PRINT 0N: ' ' 

3970 FOR 1=1 TO 24 

3980 PRINT 0N:I;*. 'jV(I) ;T(I,11;T(I,2) 

3990 NETT I 

3991 "PRINT BN:' ' 

3992 PRINT §N: "CONDITION * FDR INVERTIN6 MATRIX IS: ';Z3 
4000 RETURN 

4010 END 


.ORiGTNrr; ™ 
DF po 
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180 REM NONLINEAR DATA ANALYSIS PROGRAM (ORTHOTROPIC NQCEL) 

110 REN OETERNINE 15 OF 20 NONLINEAR COEF. BY LEAST SQUARES FIT 


120 REN USING 24 SAC MEASUREMENTS 
130 REM WRITTEN BY PROSSER 11-13-86 
140 REM 

150 REM VARIABLES 

160 REM LINEAR ELASTIC STIFFNESSES K(I> 

170 REM “ ' COMPLIANCES S(I) 

180 REM SAC MEASUREMENTS H(I,1) 

190 REM UNCERTAINTIES IN SAC'S HU, 2) 

200 REM NONLINEAR COEF' S CU,1) 

210 REM MAX. UNCERTAINTY IN NONLINEAR CQEF'S CU,2) 

220 REM PROB. " ... CCI,3) 

230 REM RECALCULATED SAC'S H(I, 3) 

240 REM MAX. UNCERTAINTY IN RECALCULATED SAC'S HU, 4) 

250 REM PROB. - ■ ■ • HU, 5) 

260 REM CQEFICIENTS OF EQ'S (24 EB. ,15 UNKNOWNS) AU,J) 

270 REM AU,J)*A(I,M) EU,J> 

280 REM INVERSE OF E(I, J> F(I,JI 

290 F.EM PRODUCT QF E AND F FIU,J) 

300 REM PRODUCT OF F AND E F2(I,J> 

310 REM CONDITION NUMBER OF E 13 

328 REM HU, 1 INCONSTANT TERMS YU) 

330 REM INTEGER COUNTERS I,J,L,H 

340 REM PERCENT DIFFERENCE BETWEEN MEASURED AND CALC SAC'S V(I) 
350 REM 1 UNCERTAINTY IN CALCULATED SAC'S T!I,2) 

1000 REM INPUT DATA 

1001 DIM K!6) ,S!6! ,H!24,5) ,A(24,15) ,CU5,3) 

1010 PA6E 


1020 PRINT ’INPUT LINEAR ELASTIC STIFFNESS MODULI (GPa) ' 

1030 PRINT ’11-1,22-2,33-3,44-4,55-5,66-6’ 

1040 FOR 1*1 TO 6 
1050 PRINT ’INPUT C’;I 
1060 INPUT KID 
1070 NEXT I 
1080 PA6E 

1090 PRINT ’INPUT LINEAR ELASTIC COMPLIANCE MODULI' (SPA! A -1’ 
1100 PRINT ’11-1,22-2,33-3,23-4,13-5,12-6’ 

1110 FOR 1=1 TO 6 
1120 PRINT “INPUT S’; I 
1138 INPUT SCI) 

1140 NEXT I 
1150 PAGE 

1160 PRINT ’INPUT SAC'S AND UNCERTAINTIES (GPa)'-l’ 

1170 PRINT ’122-1,211-2,322-3,311-4,121-5,212-6,312-7,321-8’ 

1180 FOR 1=1 TO 8 

1190 PRINT ’INPUT SAC ’;I 

1200 INPUT HI 1,1) 

1210 PRINT "INPUT UNCERTAINTY’ 

1220 INPUT HU, 2) 

1230 NEXT I 
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1248 PASE 

1258 PRINT 'SAC'S CONT'D." 

1268 PRINT '123-9,213-18,323-11,313-12,131-13,232-14,132-15,231-16* 

1278 FOR 1*9 TO 16 

1288 PRINT 'INPUT SAC ';I 

1298 INPUT H(I,1) 

1308 PRINT 'INPUT UNCERTAINTY* 

1310 INPUT H ( I ,2) 

1320 NEXT I 
1330 PAGE 

1340 PRINT 'SAC'S CON'T." 

1350 PRINT '011-17,022-18,812-19,021-20,813-21,823-22,031-23,032-24' 
1360 FOR 1=17 TO 24 
1370 PRINT 'INPUT SAC '; I 
1380 INPUT H ( 1 , 1 ) 

1390 PRINT 'INPUT UNCERTAINTY' 

1400 INPUT H ( I ,2) 

1410 NEXT I 

1420 REN CALCULATE COEF'S OF EQ.S 
1430 A=0 

1440 AU,4)=SU1/(2*K(2)> 

1441 A(1,5)=S(6)/(2*K<2)> 

1442 A(1,6)=S(5!/(2*K(21) 

1443 A(2,1)=S(6)/(2*K(1)) 

1444 A(2,2)*S(2)/(2«K(D) 

1445 A(2,3)=S(4) / (2*K ( 1 ) ) 

1446 A(3,4)=S(5I/(2*K(2! ) 

1447 A!3,5)=S(4)/(2»K(2) ) 

1448 A(3,6!=S(3)/(2*K(2)) 

1449 A (4 , l >=S<5) / (2*K(i) ) 

1450 A(4,2)=S(4I/(2*K(11) 

1451 A(4,3)=S(3)/(2*K(1)) 

1452 A(5,13)=S(1)/(2*K(6>) 

1453 A(5,14)=S(6)/(2*K(6)I 

1454 A(5,15)=S(5)/(2*K(6)1 

1455 A(6,13)=S(6)/(2*K(6>) 

1456 A(6,14I=S<2)/(2*K(6>) 

1457 A(6,15)=S(4)/!2*K(6)) 

1458 A(7,13I=S(5)/(2*K(6M 

1459 A(7,14)*S(4)/(2»K(6)I 

1460 A(7,15)=S(3)/(2*K(6)1 

1461 A(8,13)*A(7,13) 

1462 A(8,14)*A(7,14) 

1463 A(8,15)*A(7,15) 

1464 A(9,7)=S(1)/(2*K(4>) 

1465 A(9,8)=S(6)/(2«K(41) 

1466 A(9,9)=S(5)/(2*K(41! 

1467 A!10,10)=S(6I/(2*K(5! ) 

1468 A(10,11)=S(2)/(2*K(5)1 

1469 A(10,12)=S(41/(2*K(51) 

1470 A(11,7)=S(5)/(2*K(4)) 

1471 A!11,8)=S(4)/(2*K(4) ) 

1472 A(ll,9)*S(3)/(2»K(4)) 
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1473 A(12,18)=S(5)/(2*K<51) 

1474 A(12,tl)=S(4)/(2*K(5)l 

1475 A(12,12)=S(3)/(2*K(5I) 

1476 A(13,1B)=S(1)/(2»K(5)1 

1477 A(l3,ll)=S(6)/(2*K(5)) 

1478 AU3,12)=S(5)/(2*K(5)) 

1479 A(14,7)=S(6)/(2*K(4) ) 

1488 A(14,8)=S(2)/(2*K(4)I 
1481 AU4,9)*S(4)/(2»K(4)1 
1498 A(15,7)=A<9,7) 

1588 AU5,8)=A(9,8) 

1518 A(15,9)=A(9,9) 

1528 A( 16, 18)=A< 18, 18) 

1538 A(l6,ll)=A(18,ll) 

1548 A( 16, 12) =A ( 18, 12) 

1558 81=S(1I+S(6)+S(5) 

1568 B2=S<6) +S(2) +S (4) 

1578 B3=S(5)+S(41+S<3) 

1588 A(17,ll=Bl/(2*K(ll) 

1598 AU7,2)=82/{2*K(1)I 
1688 A( 17,3) =B3/ (2*K (11) 

1618 A ( 18 , 4) =B1/ (2*K (2) ) 

1628 A(18,5)=B2/(2»K(2)) 

1638 A ( 18 ,6) =B3/ (2*K ( 2) ) 

1648 A(19,13)»Bl/(2*K(6») 

1658 A<19,14)=82/(2*K(6)> 

1668 A ( 19, 15) =B3/ (2*K (6) ) 

1678 A(2B,13)=A(19,13) 

1688 A(28,14)=A(19,14) 

1698 A(28,15)=A(19,15) 

1788 A(21,1B)=81/(2*K(51) 

1718 A(21,11)=B2/(2*K(51) 

1728 A121,12)=B3/(2*K(51) 

1738 A122, 71=81/ (2»K(4) ) 

1748 A(22,8)=B2/(2»K(4)1 
1758 A(22, 91=83/ (2*K(4) 1 
1768 A(23, 18) =A (21 , 18) 

1778 A(23,1U=A(21,111 
1788 A(23,12)=A(2l,l21 
1798 A(24,7)=A(22,7) 

1888 A(24,8)=A<22,8) 

1818 A(24,9)=A(22,9) 

1888 REN AOD CONSTANT TERRS TO SAC'S 
1898 DIN Y (24) 

1988 Y(1)=-(H(1,1)+S(6) 1 
1918 Y(2)=-(H (2, 1 ) +S1.6) ) 

1928 Y (3) =-(H (3, 1 ) +S (4) ) 

1938 Y(4>=-(H(4,1)+S(5) 1 
1948 Y(5)=-(H(5,1)+S(1) 1 
1958 Y(6)=-(H(6,1)+S(2)) 

1968 Y(7)=-(H(7, 1 ) +S (4) ) 

1978 Y(8)=-(H(8,1)+S(5) 1 
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1988 Yt9)*-(H(9,1)+S(5H 
1998 Y(1B)*-<H(18,1)+S(4)) 

2888 Y(tl)=-(H(11,1)+S(3) ) 

2018 Y(12>=-CH(12 I 1)+S<3>I 
2020 Y(13)=-(H(13,1)+SU)) 

2838 Y(14)=-(H(14,1)+S(2) ) 

2848 Y(15)=-(H(15,1)+S(6) ) 

2858 Y(16)=-(H<16,1)+S16>) 

2868 Y(17I®-(H(17,1)+B1+1/(2*K(1)) ) 

2078 Y(1B) S -(HU8,1)+B2+1/(2*K(2)1) 

2080 Y(19I=-(H(19,1)+B2+1/(2*K(6))) 

2098 Y(28)=-(H(20,1)+B1+1/(2*K(6))) 

2188 Y(21)*-<H(21,1)+B3+1/<2*K(5)>) 

2110 Y(22)=-(H(22,1)+B3+1/(2»K(4))I 
2120 Y (23) =- (H (23 , l ) +B1+ 1 / (2*K (5) ) ) 

2130 Y(24)=-(H(24,1)+B2+1/(2»K(4)) ) 

2140 REM COMPUTE A(I,N)«A(I,J) 

2150 DIM E ( 15, 15) 

2160 FOR H=1 TO 15 

2170 FOR J=l TO 15 

2180 E(M,J)=0 

2190 FOR 1=1 TO 24 

2200 E(M,J)=E(N,J)+A(I,M)*A(I,J) 

2210 NEXT I 
2220 NEXT J 
2230 NEXT M 

2240 REM COMPUTE INVERSE OF E(M,J) 

2250 DIM F( 15, 15) 

2260 F=INV(E) 

2270 REM COMPUTE CONDITION NUMBER 

2280 Zl=-1 .BE+188 

2290 Z2=-l.0E+100 

2300 FOR 1=1 TO 15 

2310 FOR J=1 TO 15 

2320 IF Ed.JXZl THEN 2340 

2330 Z1=E(I,J) 

2340 IF F(I , JXZ2 THEN 2360 

2350 Z2=F(I,J1 

2360 NEXT J 

2370 NEXT I 

2380 Z3=A8S(Z2*Z1) 

2440 REN COMPUTE NONLINEAR COEF'S 

2450 FOR 1=1 TO 15 

2468 C(I,1)=0 

2470 C(I,2)=0 

2480 C ( I ,3) =0 

2490 FOR J=l TO 24 

2500 FOR 1=1 TO 15 

2510 C(I,l)=C(I,l)+F(L,I)*A(J,L)tY(J) 

2520 C(I,2)=C(I,2I+ABS(F(L,I)«A(J,L))«H(J,2) 
2530 NEXT L 
2540 NEXT J 

2550 C(I,3)=C(I,2)/23 A 0.5 
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Z56e NEXT I 

2571 REN RECALCULATE SAC VALUES USIN6 DERIVED VALUES 

2588 H<1,3)=-(2«K(2)*S(6)+S<1)*C(4,1)+S(6)»C(5,1)+S(5)*C(6,1)I 

2581 H(i t 3)*H(l,3)/(2*K(2>» 

2582 H(2,3)=-(2«K(l)*S(6)+S(6)*C(l,l)+S(2)*C(2,t)+S(4)*C(3,l)) 

2583 H(2,3)=H(2,3)/(2*K(1)) 

2584 H(3,3)=-(2*K(2)*S(4)+S(5)*C(4,l)+S(4)*C(5,l)+S(3)tC(6,l)) 

2585 H(3,3)=H(3,3)/(2«K(2)) 

2586 H(4,3)=-(2*K(l)»S(5)+S(5)*Cll,l)+S(4)tC(2,l)+S(3)*C(3,l>) 

25B7 H(4,3)=H(4,3)/(2«K(1)) 

2588 H(5,3) s -(2*K(6l*S(l)+S(l)*Ctl3,l)+S(6)*C(14,l)+S(5)*C(15,l) ) 

2589 H(5,3)=H(5,3)/(2*K<6)) 

2591 H(6,3)=-(2*K(6)*S(2)+S(6)»C(13,l)*S(2)«C(l4,t)+S(4)»C(15,ll) 

2591 H(6,3)=H(6,3I/(2*K(6)) 

2592 H(7,3)=-(2«K(6)*S(4)+S(5)*C(13,l)+S(4)*C(l4,t)+S(3)*C(15,l)) 

2593 H(7,3)=H(7,3I/(2«K<6)) 

2594 H(0,3)=-(2*K<6I*S(5)+S<5)*C<13,1>+S(4)«C(14,1)+S(3)*C(15,1)) 

2595 H(8,3)*H(8,3)/(2#K(6)) 

2596 H(9,3)=-(2*K(4)#S(5)*S(l)*C(7,l)+S(6)*C(8,l>+S(5>»C(9,l>) 

2597 H(9,3)=H(9,3)/(2«(4)) 

2598 H(18,3)=-(2»K(5)*S(4)+Sl6)*C(18,l)+S(2)*C(ll,l)+S(4)tC(12,l)) 

2599 H(18,3)=H(1B,3)/(2*K(5)) 

2688 H(ll,3) s -(2»K(4)*S(3)+S(5)*C(7,l)+S(4)*C(8,l)+S(3)*C(9,l) ) 

2601 H(11,3) S H(U,3)/(2*K(4) ) 

2602 H(12,3) s -(2*K(5)*S(3)+S(5)*C(l0,l)+S(4)*C(ll,l)*S(3)*C(!2,l!) 

2683 H(12,3)=H(12,3) /(2*K(5) ) 

2684 H(13,3)*-(2*K(5)*S(lH-S(l)«C(10,l)+S(6)«C;il,i)+S!5!*C(12,l>) 
2605 H(13,3)=H(13,3>/(2«K(5)) 

2686 H(t4,3)=-(2*K(4)*S(2)+S(6)*C<7,l)+S!2>»C(8,l)+S!4)*C(9,l)) 

2687 H(14,3)=H(14,3)/(Z«K(4!) 

2688 H(l5 I 3)=-(2*K(4)*Si6)+S(l)*C(7 1 l)+S(6)*C(8 l l)+S(5HC(9,l)) 

2609 H(15,3)=H(15,3)/(2*K(4) ) 

2610 H(16,3I=-(2*K(5)*S(6)+S(6)*C(10,1)+S(2)*C(11,1I+S(4>*C(12,1)) 

2611 H(16,3)=H(16,3) / (2«K(5) ) 

2620 H(l7 f 3)=-(l+2*Bl*K(l)+Bl«C(l,l)+B2*C<2,l)+B3*Ci3,l!) 

2630 H(17,3)=H(17,3)/(2*K(1)) 

2640 H(18,3)=-(t+2»K(2)*B2+Bl*C(4,l)+B2*C(5,l)+B3*C(6,l)) 

2650 H(IB,3)=H(18,3)/(2*K(2)) 

2660 H ( 19,3) =- ( 1+2*K (6) *B2+B1*C (13,1) +B2*C (14,1) +B3*C (15,1) ) 

2670 H(19,3)=H(19,3)/(2*K(6)) 

2680 H (20, 3) =- ( 1+2*K (6) *B1+B1*C( 13,1) +B2»C( 14, 1) +B3*C (15,1) ) 

2690 H(28,3)=H(20,3) / (2*K(6) ) 

2708 H (21 , 3) =- ( 1+2*K (5) *B3+B 1*C (18,1) +B2*C (11,1) +B3*C (12,1)) 

2710 H(21,3)=H(21,3)/(2*K(5) ) 

2728 H(22,3)=-(1+2*K(4)*83+B1*C(7,1I+B2*C(B,1H-B3«C(9,1)) 

2730 H(22,3) =H (22,3) / (2*K (4) ) 

2748 H(23,3)=-(l+2*K(5)*Bl+Bl*C(10,l)+B2*C(ll,l)+B3»C(12,l) ) 

2750 H(23 ,3) =H (23 ,3) / (2*K (5) ) 

2768 H (24,3)=- ( l+2*K (4) *B2+Bl»C(7, 1 ) +B2*C (8, 1 )+B3*C(9 ,1)1 
2770 H(24,3)=H(24,3)/(2*K(4)) 

2810 DIR 0(15,2) 

2820 FOR 1=1 TO 2 
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2830 FOR J=l TO 15 

2840 W(J , I )-ABS <C ( J , 1+1 ) /C ( J , 1 i i 

2950 NEXT J • 

2960 NEXT I 
2970 FOR 1=1 TO 2 

2880 H(l,I+3)=iN(4,I)+N(5,Il+N(6,l) !*H(l,3) 

2881 H(2,I+3)=!N(l,I)+N(2, I)*N(3,I) i*H(2,3) 

2892 H(3,I+3)=!«(4,I>+N(5,I)+N(6,I)!*H(3,3) 

2883 H(4, 1+3) = (N<1, 11*11(2, 11+0(3, 1) !*H(4, 3) 

2884 H<5,I+3)=(W(13,I)+N(14,!)+N<15,I!)*H(5,3) 
2985 H(6, I+3)=(M(13, I1+IM14, 11+0(15, I) )«H(6, 31 

2886 H(7,I+3) = !N(13,I!+«(l4,I!+ii(15,I!i*H(7,3) 

2887 H(8,I+3) = (N(13,I!+N(14,n+N(t5,I!>*H(8,3) 

2888 H(9,I+3) = (i4(7,i>+0<8,l)+N(9,Ill*H(9,3) 

2889 H(10,I+31=(N(10,I!+N(ll,I)+H(l2,I) I*H(10,3) 

2890 H(ll,I+3! = (W(7,I)+«(8,I>+IK9,I)}*H(ll,3> 

2891 HU2,i+3) = (H<18, 11+0(11, 11+0(12, 1))*H(12, 3) 

2892 H( 13, 1+3! =(0(10,11+0(11, 1) +0(12, 1) I *8(13,3) 

2893 H(14,I+3)=(0(7, 11+0(8, 11+0(9, I))*H(14, 3) 
2994 H(15,I+3)=(0(7,I)+0(8,I)+0(9,I!1*H(15,3) 

2895 H(16,I+3)=(0(10,I)+0(11,I)+0(12,I) i*H(16,3) 

2896 H(!7,I+3!=(0(l,l!+0(2,!)+0(3,I) >+H(!7,3) 

2897 H(19,I+3S=(0(4,I!+0(5,I)+0(6,I! )+H(18,3) 

2898 H( 19,1+3)= (0(13, I) +0(14,11+0 (15, I) )»H( 19, 31 

2899 H (20, 1+31= -0(13, I! +0(14, 11+0(15, I! !*H!20,3) 

2900 H(21,!+3!=(0(!0,[)+0(l 1,11+0(12, I! !*H(21,3) 

2901 H!22,I+3)=(0(7,I!+0(3,I)+0(9,!))*H(22,3) 

2902 H(23,!+3)=(0(18,n+0(ll,I)+0(12,I) !*H(23,3! 

2903 H!24,I+3;=(0(7,I)+0(B,!>+0(9,!))+H(24,3) 
3860 NEXT I 

3070 FOR 1=4 TO 5 
3080 FOR 3=1 TO 24 
3090 H(J,I1=ABS!H(J,I1) 

3180 NEXT J 
3118 NEXT ! 

3128 DIN V (24) 

3133 FOR 1=1 TO 24 

3140 V ! I ! = A8S ( ( H ! 1 , 1 ) -H ( 1 , 3 ) ) / H( 1 , 1 ! ) 

3158 v’(Ii=V<I>*188 

3160 NEXT I 

3178 DIN T (24,21 

3180 FOR 1=1 TO 2 

3190 FOR J=1 TO 24 

3280 T(J,I)=H(3,I+3i/H(J, 31*108 

3201 T!J,I1=ABS!T(J,I1) 

3218 NEXT 3 
3223 NEXT ! 

3230 REN PRINT OUT DATA 

3240 REN 

3250 N=32 

3260 PAGE 

3270 80SUB 3380 

3280 PRINT "PRINTED COPY? YES =1 NO =2 " 


ORIGINAL PAGE IS 
OF POOR QUALITY 
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3300 IF Z5=l THEN 3340 
3310 IF- 15=2 THEN 4010 
3320 PRINT 'BAD CHOICE' 

3330 GO TO 3280 

3340 PRINT 'INPUT PRINTER I' 

3350 INPUT N 
33i0 30SUB 3380 
3370 GO TO 4010 

3380 PRINT AN: 'NONLINEAR DATA ANALYSIS FOR ORTHOTROPIC MODEL 
3390 PRINT «N:' " 

3400 PRINT AN:' INPUT DATA' 

3410 PRINT AN:' ' 

3420 PRINT AN: "LINEAR ELASTIC MODULI" 

3430 PRINT AN:* ' 

3440 PRINT AN:"Cli 
3450 PRINT AN:'C22 
3460 PRINT AN:*C33 

3470 PRINT AN: "C44 = 

3471 PRINT AN:"C55 = ";K(5); 

3472 PRINT AN:'C66 = “;K(6); 

3480 PRINT AN:" " 

3490 PRINT AN: 'VELOCITY DERIVATIVE DATA iGPaJM" 

3500. PRINT AN:' " 

H122 = " ; H ( 1 , l ! ; " +/- ";H!l,2i 
H211 = ’;H(2,ii:' +/- ";H(2,2) 

H322 = *;H(3,1);* ♦/- ";H(3,2: 



!3Pa! 

511 = 


!SPa) A -t 

;K(2)j * 

IGPa) 

S22 = 

";S(2I;" 

(GPa) A_ l 

;K(3) :* 

< GPa) 

S33 = 

” ; S ( 3 > ; * 

(SPa) '-I 

;K <4) ; " 

(SPa) 

323 = 

-;S(4);' 

(GPa) A -l 

;K(5);" 

(GPa) 

513 = 

■;S(5);’ 

(SPa) ‘-1 

; K (A) ; " 

(GPa) 

S12 = 

" ; S (6 ) ; “ 

(GPa)'-l 


3310 PRINT AN: '1. 
3520 PRINT AN: "2. 
3530 PRINT AN: "3. 
3340 PRINT AN: '4. 
3550 PRINT AN:"5. 
3560 PRINT AN: "6. 
3570 PRINT AN: '7. 
3580 PRINT AN: "8. 
3590 PRINT AN: "9. 


H311 


*; H<4, 1 ) ; " +/- , jH(4 1 21 


H121 = ' ;H (5, 13 ; ' ♦/- ";H(5,2) 
H212 = ' ; H ( 6 , 1 ) ; " ♦/- ";H!6,2; 
H312 = •;H(7,Uj* ♦/- “;HI7,2> 
H321 = *;H(8,1);' +/- ";H!S,2) 
H123 = *;H(9,1!;' ♦/- ”;H(9,2) 


3600 PRINT 9N:*10. H21I = "jH'10,1!;" +/- ";H<10,2 


3610 PRINT AN:‘U. 

H323 

3 

•}H(11,1>5" 

4* / - 

";H!11,2.) 

3620 PRINT AN: *12. 

H313 

= 

";H(12,D ; ” 

+ / - 

";H(12,2i 

3630 PRINT AN: *13. 

H131 

= 

*;H!13,U;" 

4 ./- 

*» . L| ,M 7 T, 
, fl •. i Jji 

3640 PRINT AN: '14. 

H232 

= 

';H(14,H 

+ /- 

*;H(14,2) 

3650 PRINT AN: "15. 

H132 

3 

*;H(15,1! 

+ /“* 

";H(15,2! 

3660 PRINT AN: "16. 

H231 

= 

*;H(16,li:“ 

+ /- 

* : H ( 1 6 , 2 ) 

3672 PRINT AN:' 17. 

mil 

= 

";H(17,1!;" 

+ /- 

*;H(17,2) 

3680 PRINT AN: "19. 

H022 

3 

*;H(18,15;" 

t/“ 

“;H(1S,2! 

3690 PRINT AN: '19. 

H012 

= 

*;H(19,U;" 

+ /- 

";H(19,2! 

3700 PRINT AN: "20. 

H021 

= 

*;H(20,l)j" 

+ /* 

' ; H (20,2) 

3710 PRINT AN: "21. 

H013 

3 

’;H(2l,i)j" 

♦ /- 

■;H(21,2) 

3728 PRINT AN: "22. 

H023 

= 

' ; H (22 , l ) ; " 

+ /- 

" ; H ( 22 ,2) 

3732 PRINT AN: ”23. 

H031 

3 

•;H(23,1);" 

+ /~ 

";H(23,2! 

3740 PRINT AN: '24. 

H232 

3 

*;H(24,l); " 

t/» 

";K(24,21 


3750 PRINT AN:' " 

3760 PRINT AN: 'THIRD ORDER ELASTIC CONSTANTS DERIVED FROM LEAST' 
3770 PRINT AN: 'SQUARES FIT AND MAX. ERROR AND PRQB. ERROR (GPa)' 
3780 PRINT AN:' " 
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ORIGINAL PAGE 13 
OF. POOR QUALITY 


3798 PRINT 0N:*81ll 
3808 PRINT 8N;*C112 
3810 PRINT SN:*C113 
3820 PRINT §N: "C122 
3830 PRINT BN: *C222 
3840 PRINT BN:*C223 

3850 PRINT BN:*C144 

3851 PRINT BN: "8244 

3852 PRINT BN: *8344 

3853 PRINT SN:"C155 

3854 PRINT BN:*C255 

3855 PRINT BN:'C355 

3856 PRINT SN:'C166 

3857 PRINT BN: *0266 
3353 PRINT BN: *0366 


■sCd.llj* +/- * 58 ( 1 , 2 ) 5 * ♦/- * 58 ( 1,31 

*;8<2,1);' ♦/- -;C(2 t 2);* ♦/- *;C(2,3) 

*;C(3, 1) ; * +/- *;8(3,2);* ♦/- *;C(3,3) 

* ; C (4, 1 ) ; * ♦/- *;C(4,2);' +/- *;C!4,3) 

* ; C (5, 1 ) ; * ♦/- ";C (5,2) ; * +/- *;C(5,3) 

•. 58 ( 6 , 1 );* ♦/- -;C(6,2);* ♦/- *;8(6,3) 

" ; C ( 7 , 1 ) ; * +/- *;8(7,2);* +/- *;C(7,31 

*;8(8,1);' ♦/- , ;C(8,2);* +/- * 58 ( 8 , 3 ) 

-;C(9,1) ; * +/- * 58 ( 9 , 2 ) 5 * +/- *;C(9,3) 

*;8(1B, 1) ; * ♦/- *;C(10,2);* +/- ";C(10,31 
*58(11,1);' ♦/- *;8(il,2);* +/- *;C(ll,3) 
*58(12,1)5* ♦/- *58(12, 2);* +/- *;8(12,3) 
*58(13,1)5* ♦/- * ;C( 13,2) ; * ♦/- *;C(13,3) 
*58(14,1)5* ♦/- *;C(14,2);* +/- *;Cil4,3) 
*;8(15,t);* ♦/- *;C(15,2>;* +/- *;C(15,3) 


3860 PRINT 8N: * 8 133 , C233 , C333 , 8 1 23 , C456 UNDETERMINED * 

3870 PRINT BN:* * 

3880 PRINT BN: 'RECALCULATED VELQ8ITY DERIVATIVES USING FITTED TOEC'S* 
3890 PRINT SN:*AND MAX. AND PROB. ERRORS IN SANE ORDER AS BEFORE* 


3900 PRINT BN:* " 


3910 FOR 1=1 TO 24 

3922 PRINT 6N:I;\ *;HU,3!;* ♦/- *;H(I,4);* ♦/- *;H(I,5) 
3930 NEXT I 


3940 PRINT BN: * " 

3950 PRINT SN: "PERCENT DIFFEREN8ES, NAX. AND PROB. UNOERTAINTIES* 
3960 PRINT BN: * * 

3970 FOR 1=1 TO' 24 


3980 PRINT BN: I;*. *;V(I);T(I,ll;T(l,2) 
3990 NEXT I 


3991 PRINT BN: * * 

3992 PRINT BN: 'CONDITION 4 FOR INVERTING MATRIX IS: *; 23 
4000 RETURN 

4010 ENO 



- 184 - 




REFERENCES 


[1] F. D. Murnaghan, "Finite Deformation of an Elastic Solid," Am. J. Math , 49 (1937) 
pp. 235-260. 

[2] F. Birch, "The Effect of Pressure Upon the Elastic Parameters of Isotropic Solids, 
According to Murnaghan’s Theory of Finite Strain," J. Appl. Phys., 9 (1938) pp. 
279-288. 

[3] M. Born, Handbuch dtr Physik XXIV/2 , (Springer, Berlin, 1933) pp. 629-724. 

[4] E. Madelung, Zeits fur Physik 11 (1910) p. 898. 

[5] R. E. Green, Jr., "Ultrasonic Investigation of Mechanical Properties," Treatise on 
Materials Science and Technology , Vol. 3 , H. Herman, ed., (Academic Press, New 
York, 1973) pp. Ml, 73. 

[6] A. A. Gedroits and V. A. Krasil’Nikov, "Finite Amplitude Elastic Waves in Solids 
and Deviations from Hooke’s Law," Soviet Physics J.E.T.P., 16 (1963) pp. 1122- 
1126. 

[7] M. A. Breazeale and D. O. Thompson, " Finite Amplitude Ultrasonic Waves in 
Aluminum," Appl. Phys . Let., 3 (1963) pp. 77-78. 

[8] F. R. Rollins, Jr., "Interaction of Ultrasonic Waves in Solid Media," Appl. Phys . 
Let., 2 (1963) pp.147-148. 

[9] F. R. Rollins, Jr., L. H. Taylor, and P. H. Todd, Jr., "Ultrasonic Study of Three- 
phonon Interactions," Phys. Rev., 136 (1964) pp. A597-A601. 

[10] D. S. Hughes and J. L. Kelly, "Second Order Elastic Deformation of Solids," Phys. 
Rev., 92 (1953) pp. 1145-1149. 

[11] R. A. Toupin and B. Bernstein, "Sound Waves in Deformed Perfectly Elastic 
Materials, Acoustoelastic Effect," J.A.S.A., 33 (1961) pp. 216-225. 

[12] K. Brugger, "Thermodynamic Definition of Higher Order Elastic Coefficients," 
Phys. Rev., 133 (1964) pp. A1611-A1612. 

[13] K. Brugger, "Determination of Third-Order Elastic Coefficients in Crystals," J . 
Appl . Phys., 36 (1965) pp. 768-773. 

[14] R. N. Thurston and K. Brugger, "Third-Order Elastic Constants and the Velocity 
of Small Amplitude Elastic Waves in Homogeneously Stressed Media," Phys. Rev., 
133 (1966) ppA1604-A1610. 

[15] D. S. Hughes, E. B. Blankenship, and R. L. Mims, ’Variation of Elastic Moduli and 
Wave Velocity with Pressure and Temperature in Plastics," J. Appl. Phys., 21 
(1950) pp. 294-297. 

[16] H. J. McSkimin, "Measurement of Elastic Constants at Low Temperatures by 
Means of Ultrasonic Waves - Data for Silicon and Germanium Single Crystals, and 
for Fused Silica," J. Appl. Phys., 24 (1953) pp. 988-997. 

[17] H. J. McSkimin, "Use of High Frequency Ultrasound for Determining the Elastic 
Moduli of Small Specimen," Presented at the National Electronics Conference, 
1956, Chicago, Illinois. 

[18] J. Williams and J. Lamb, "On the Measurement of Ultrasonic Velocity in Solids," 
J.A.S.A. , 30 (1958) pp. 308-313. 

[19] H. J. McSkimin, "Pulse Superposition Method for Measuring Ultrasonic Wave 
Velocities in Solids," J.A.S.A ., 33 (1961) pp. 12-16. 


i 


c -7> 



- 185 - 


[20] J. E. May, Jr., L R. E. Nat . Conv . Rec., 6 (1958) pp. 134-142. 

[21] E. P. Papadakis, "Ultrasonic Phase Velocity by the Pulse-Echo-Overlap Method 
Incorporating Diffraction Phase Corrections," J.A.S.A. , 42 (1967) pp. 1045-1051. 

[22] E. P. Papadakis, "Absolute Accuracy of the Pulse-Echo-Overlap Method and the 
Pulse-Superposition Method for Ultrasonic Velocity," J.A.S.A., 52 (1972) pp. 843- 
846. 

[23] J. S. Heyman and Francis D. Stone, "Pulsed Phase Locked Loop Strain Monitor," 
NASA Tech Briefs , 6, (1981) p. 68. 

[24] J. S. Heyman and E. J. Chern, "Ultrasonic Measurement of Axial Stress," J. of 
Testing and Evaluation, 10 (1982) pp. 202-211. 

[25] N. G. Einspruch and R. J. Manning, "Third-Order Elastic Moduli of Anisotropic 
Solids," J . Appl. Phys 35 (1964) pp.560-567. 

[26] R. W. Benson and V. J. Raelson, "Acoustoelasticity," Product Eng., 30 (1959) pp. 
56-59. 

[27] N. N. Hsu, "Acoustical Birefringence and the Use of Ultrasonic Waves for Experi- 
mental Stress Analysis," Exp. Mech., 14 (1974) pp. 169-176. 

[28] Y. Pao, W. Sachse, and H. Fukuoka, "Acoustoelasticity and Ultrasonic Measure- 
ments of Residual Stresses," Physical Acoustics Vol. XVII., W. P. Mason and R. N. 
Thurston, ed., (Academic Press, New York, 1984) pp. 62-140. 

[29] A. V. Clark, R. B. Mignona, and R. J. Sanford, "Acousto-elastic Measurement of 
Stress and Stress Intensity Factors Around Crack Tips," Ultrasonics , 21 (1983) 
pp. 57-64. 

[30] G. S. Kino, J. B. Junter, G. C. Johnson, A. R. Selfridge, D. M. Barnett, G. Her- 
mann, and C. R. Steele, "Acoustoelastic Imaging of Stress Fields," J. Appl. Phys., 
50 (1979) pp. 2607-2613. 

[31] J. S. Heyman, S. G. Allison, and K. Salama, "Influence of Carbon Content on 
Higher-Order Ultrasonic Properties in Steels," IEEE Ultrasonics Symposium, (1983) 
pp. 991-994. 

[32] J. S. Heyman, and E. J. Chern, "Characterization of Heat Treatment in Aluminum 
Based on Ultrasonic Determination of the Second and Third Order Elastic Con- 
stants," IEEE Ultrasonics Symposium, (1981) pp. 936-939. 

[33] H. Singh and A. W. Nolle, "Pressure Dependence of the Viscoelastic Behavior of 
Polyisobutylene," J. Appl. Phys., 30 (1959) pp. 337-342. 

[34] J. R. Assay, D. L. Lamberson, and A. H. Guenther, "Pressure and Temperature 
Dependence of the Acoustic Velocities in Polymethylmethacrylate," J. Appl. Phys., 
40 (1969) pp. 1768-1783. 

[35] J. R. Assay, D. L. Lamberson, and A. H. Guenther, "Ultrasonic Technique for 
Determining Sound Velocity Changes in High-Loss Materials," J.A.S.A., 45 (1969) 
pp. 566-571. 

[36] D. L. Lamberson, Ph.D. dissertation, Air Force Institute of Technology, 1969. 

[37] L. K. Zarembo and V. V. Shklovskaya-Kordi, "Generation of Harmonics During 
Propagation of Ultrasonic Longitudinal and Shear Waves in Solids," Soviet Physics 
- Solid State , 12 (1971) pp. 2960-2961. 

[38] Y. Wada, A. Itani, T. Nishi, and S Nagai, "Gruneisen Constant and Thermal Pro- 
perties of Crystalline and Glassy Polymers," J. Poly. Sci ., 7 (1969) pp. 201-208. 

[39] D. L. Lamberson, J. R. Assay, and A. H. Guenther, "Equation of State of Polys- 
tyrene and Polymethylmethacrylate from Ultrasonic Measurements at Moderate 
Pressures," J. Appl. Phys., 43 (1972) pp. 976-985. 



- 186 - 


[40] R. E. Hankey and D. E. Schuele, "Third Order Elastic Coefficients of A1 2 0 3 
J.A.S.A., 48 (1970) pp. 190-202. 

[41] E. H. Love, On the Mathematical Theory of Elasticity, (Cambridge University Press, 
New York, 1944). 

[42] K. Brugger, "Pure Modes for Elastic Waves in Crystals," J. Appl. Phys., 36 (1965) 
pp. 759-768. 

[43] W. P. Mason, Physical Acoustics and the Properties of Solids, (Nostrand Co., Inc., 
Princeton, 1958) pp. 95-97. 

[44] L. Y. Tu, J. N. Brennan, and J. A. Sauer, "Dispersion of Ultrasonic Pulse Velocity 
in Cylindrical Rods," J.A.S.A., 27 (1955) pp. 550-555. 

[45] R. D. ICriz and W. W. Stinchcomb, "Elastic Moduli of Transversely Isotropic Gra- 
phite Fibers and their Composites," Experimental Mechanics, 19 (1979) pp. 41-49. 

[46] D. Hull, "An Introduction to Composite Materials," Cambridge Solid State Science 
Series, R. W. Cahn et al., ed., (Cambridge University Press, New York, 1981) p. 7. 

[47] H. M. Ledbetter and D. T. Read, "Orthorhombic Elastic Constants of a NbTi/Cu 
Composite Superconductor," J.A.P., 48 (1977) pp. 1874-1879. 

[48] D. P. Garber, "Tensile Stress-Strain Behavior of Graphite/Epoxy Laminates,” 
NASA Contractor Report 8592 Contract NAS1-16000, (1982) p. 21. 

[49] J. H. Cantrell, Jr., "Anharmonic Properties of Solids from Measurements of the 
Stress Acoustic Constant," J. of Testing and Evaluation, 10 (1982) pp. 223-229. 


Standard Bibliographic Page 


1. Report No. 

NASA CR-4100 


2. Government Accession No. 


3. Recipient’s Catalog No. 


j 4. Title and Subtitle 5. Report Date 

Ultrasonic Characterization of the Nonlinear Elastic October 1987 
Properties of Unidirectional Graphite/Epoxy Composites 6 . Performing Organization Code 


7. Author(s) 


William H. Prosser 


8. Performing Organization Report No. 


10. Work Unit No. 


11. Contract or Grant No. 


9. Performing Organization Name and Address 

Center for NDE 506-43-11-03 

102 Maryland Hall 11. Contract or Grant No. 

pe Johns Hopkins University NGT 21-001-802 

Baltimore, MD 21218 — _ — — — 


12. Sponsoring Agency Name and Address 

National Aeronautics and Space Administration 
Langley Research Center 
Hampton, VA 23665-5225 


15. Supplementary Notes _ n „ 

Langley Technical Monitor: Joseph S. Heyman 

Final Report 

The information presented in this report was offered as an essay in conformity with 
the requirements for the Master of Science in Engineering degree, The Johns Hopkins 
University, Baltimore, Maryland, March 1987. 


NGT 21-001-802 

13. Type of Report and Period Covered 

Contractor Report 


14. Sponsoring Agency Code 


i 


16. Abstract 

Because of their strength, stiffness, and light weight, fiber reinforced composite 
materials offer great potential for applications particularly in the aerospace 
industry. Early research has proven the usefulness of these materials, but the neei 
to quantitatively characterize important material properties and develop applicable 
nondestructive evaluation techniques remains. Linear and nonlinear elastic 
properties are important physical properties in conventional materials in this 
respect. Nonlinear properties are important in the nondestructive determination of 
applied and residual stress as well as measuring the interatomic bonding forces in 
crystalline solids. Also, several investigations have established a possible 
relationship between nonlinear elastic properties and ultimate strength in aluminum 
and carbon steel. 

This work presents the theoretical treatment of linear and nonlinear elasticity in . 
unidirectionally fiber reinforced composite as well as measurements for a 
unidirectional graphite/epoxy composite (T300/5208). Linear elastic properties wer< 
measured by both ultrasonic and strain gauge measurements. The nonlinear properties 
were determined by measuring changes in ultrasonic "natural" phase velocity with a 
pulsed phase locked loop interferometer as a function of stress and temperature. 
These measurements provide the basis for further investigations into the 
relationship between nonlinear elastic properties and other important properties 
such as strength and fiber-matrix interfacial strength in graphite/epoxy composites 
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